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Abstract 

A  key  application  of  Laser  Detection  and  Ranging  (LADAR)  systems  is  measure¬ 
ment  of  range  to  a  target.  Many  modern  LADAR  systems  are  capable  of  transmitting 
laser  pulses  that  are  less  than  a  few  nanoseconds  in  duration.  These  short-duration 
pulses  provide  excellent  range  precision.  However,  randomness  in  the  detected  laser 
signals  places  limits  on  the  precision. 

The  goal  of  this  dissertation  is  to  quantify  the  range  precision  limits  of  LADAR 
systems.  The  randomness  in  the  time  between  photon  arrivals,  which  is  called  shot 
noise,  is  discussed  in  depth.  System-dependent  noise  sources  such  as  dark  current  and 
detector  gain  variation  are  considered.  The  effect  of  scene- dependent  parameters  in¬ 
cluding  background  light,  target  obscuration,  and  target  orientation  is  also  discussed. 
Finally,  noise  mitigation  strategies  such  as  pulse  averaging  and  gain  equalization  are 
described  and  tested  on  simulated  and  real  LADAR  data. 
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OE  LADAR  Systems 

I.  Introduction 

Laser  detection  and  ranging  or  LADAR  systems  have  developed  rapidly  in  recent 
years.  LADAR  is  analogous  to  radar,  but  is  done  at  optical  wavelengths  rather  than 
radio  frequencies.  The  use  of  optical  wavelengths,  which  are  roughly  one-thousandth 
the  length  of  the  shortest  radar  wavelengths,  yields  many  benefits  and  challenges. 

LADAR  systems  are  used  for  many  remote  sensing  applications.  They  are 
used  to  measure  wind  speeds  in  clear  air  [19].  LADAR  systems  with  two  or  more 
wavelengths  are  used  to  detect  and  measure  concentrations  of  certain  gasses  in  the 
atmosphere  [28].  Some  LADARs  interfere  received  laser  light  with  a  local  laser  in 
order  to  measure  Doppler  shift  [8].  These  coherent  LADARs  are  sensitive  enough  to 
measure  surface  vibrations  on  remote  objects  [17].  However,  a  primary  application  of 
LADAR  has  always  been  the  measurement  of  range  to  a  target. 

Different  LADAR  technologies  are  used  to  measure  target  range  with  a  laser. 
Amplitude  and  frequency  modulation  of  continuous-wave  lasers  has  been  used  to 
measure  range  [44].  But  the  most  commonly-used  technique  is  to  transmit  a  short- 
duration  laser  pulse  and  measure  its  time  of  flight.  Multiple  range  measurements 
taken  along  different  lines  of  sight  can  be  combined  to  create  a  three-dimensional 
image  of  a  scene.  These  three-dimensional  images  are  valuable  in  applications  such 
as  mapping  and  target  recognition. 

The  precision  of  range  measurements  made  with  a  LADAR  system  is  limited  by 
noise.  The  fundamental  source  of  noise  in  a  laser  signal  is  the  random  time  between 
photon  arrivals.  This  phenomenon,  called  shot-noise,  can  be  modeled  as  a  Poisson 
random  process  [21].  Other  noise  sources  in  LADAR  systems  can  also  be  modeled  as 
Poisson  processes  [31]. 
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Noise  in  LADAR  systems  is  not  always  modeled  using  the  Poisson  distribution. 
The  noise  in  systems  with  high  levels  of  electronic  noise  can  be  accurately  character¬ 
ized  with  a  Gaussian  distribution  [22] .  The  various  noise  models  (Poisson  or  Gaussian) 
can  be  used  to  derive  limits  on  the  precision  of  range  measurements. 

This  dissertation  examines  several  different  signal  and  noise  models  for  LADAR 
systems  transmitting  short  pulses  of  laser  light.  Limits  on  range  precision  are  derived 
for  different  scenarios  such  as  calculating  range  to  a  target  that  is  partially  obscured 
or  is  tilted  at  an  oblique  angle  when  the  LADAR  observes  it.  Electronic  effects  in 
systems  can  distort  the  data  recorded  by  a  LADAR.  This  distortion  can  affect  range 
measurements.  Thus,  methods  of  correction  are  developed  and  tested.  In  all  cases,  the 
goal  of  this  dissertation  is  to  characterize  the  precision  of  LADAR  range  measurements 
and  to  develop  new  methods  that  improve  the  process  of  range  estimation. 

1.1  LADAR  Systems 

1.1.1  LADAR  System  Diagram.  All  LADAR  systems  include  a  laser  trans¬ 
mitter  and  a  receiver  designed  to  detect  laser  light  reflected  from  a  target.  A  LADAR 
system  block  diagram  is  shown  in  Fig.  1.1.  The  transmitter  emits  laser  light  toward 
the  target.  In  an  imaging  system  the  target  is  the  solid  object  in  the  scene.  However, 
some  LADAR  systems  are  designed  to  observe  laser  light  scattered  from  aerosols.  A 
fraction  of  the  light  hitting  the  target  is  reflected  into  the  LADAR’s  receiver.  That 
light  is  detected  and  processed  by  the  system. 

A  wide  variety  of  lasers  have  been  used  as  LADAR  transmitters.  Wavelengths 
between  the  ultraviolet  and  the  long  wave  infrared  have  been  employed.  However, 
it  is  most  common  to  use  lasers  that  operate  in  the  near  infrared  when  building 
a  LADAR  system.  The  near  infrared  band  of  wavelengths  offers  reasonably  good 
eye-safety.  Because  most  optical  telecommunication  equipment  operates  in  the  near 
infrared,  optical  components  designed  to  operate  at  those  wavelengths  are  plentiful 
and  relatively  inexpensive. 
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Figure  1.1:  Block  diagram  used  to  model  a  LADAR  system. 


LADAR  systems  also  need  receivers  to  detect  reflected  laser  light.  The  receiver 
consists  of  a  lens  that  focuses  reflected  laser  light  onto  a  photodetector.  The  photode¬ 
tector  converts  the  received  light  into  an  electrical  signal.  That  signal  is  processed  in 
order  to  measure  information  about  the  target  such  as  range,  reflectivity,  and  velocity. 

Some  LADAR  systems  use  the  same  aperture  to  transmit  and  receive  laser  light. 
Such  systems  are  called  monostatic.  Systems  with  separate  transmit  and  receive 
apertures  are  referred  to  as  bistatic.  Bistatic  systems  usually  have  transmit  and 
receive  apertures  that  are  located  very  near  each  other. 

1.1.2  Pulsed  LADAR  Systems.  Laser  detection  and  ranging  systems  can  be 
used  to  measure  range  to  targets.  In  order  to  measure  range,  typical  LADAR  systems 
transmit  short  pulses  of  laser  light.  There  are  alternative  methods  of  measuring 
target  range  such  as  with  a  chirped  waveform  and  coherent  processing.  But  the  most 
common  method  of  range  estimation  involves  measuring  the  time  of  flight  for  a  short 
duration  laser  pulse. 

Typical  LADAR  systems  transmit  short  laser  pules  using  a  process  called  Q- 
switching.  The  Q-switching  process  allows  relatively  large  amounts  of  energy  to  build 
within  the  laser  cavity  before  being  released  over  a  very  short  period  of  time.  The 
energy  is  released  from  the  laser  cavity  using  a  device  such  as  a  rotating  mirror  or 
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a  Pockels  cell  crystal  [38].  Laser  pulses  created  by  Q-switched  systems  can  have 
durations  measured  in  nanoseconds. 

Given  a  fixed  pulse  energy,  shorter  pulses  provide  better  range  precision.  How¬ 
ever,  shorter  pulses  require  optical  components  and  photodetectors  that  can  tolerate 
higher  peak  powers.  Shorter  pulses  also  demand  faster  sampling.  Generation  of 
laser  pulses  is  often  accomplished  using  a  component  within  the  laser  cavity  called 
a  Q-switch.  Shorter  pulses  necessitate  the  use  of  more  sophisticated  and  expensive 
Q-switches.  These  engineering  issues  place  practical  limits  on  how  short  a  laser  pulse 
can  be  made.  Gurrent  LADAR  systems  typically  transmit  pulses  that  are  a  few 
nanoseconds  long  and  sample  those  pulses  at  a  rate  around  one  Gigahertz. 

The  short  pulses  created  by  the  Q-switching  process  yield  excellent  range  preci¬ 
sion.  LADAR  range  precision  of  less  than  3  cm  has  been  claimed  [2] .  However,  short 
pulses  place  demands  on  the  detection  and  digitization  hardware  in  LADAR  systems 
since  nanosecond-class  pulses  have  bandwidths  on  the  order  of  a  Gigahertz. 

1.1.3  Imaging  LADAR  Systems.  Development  of  LADAR  imaging  systems 
is  advancing  rapidly.  Early  LADAR  imaging  systems  used  few  detectors  and  usu¬ 
ally  relied  on  mechanical  scanning  systems  to  produce  images  with  large  numbers 
of  pixels.  Newer  systems  use  arrays  of  photodetectors  and  can  produce  images  with 
large  numbers  of  pixels  without  the  use  of  mechanical  scanning.  Imaging  LADARs 
with  more  than  ten  thousand  pixels  are  commercially  available.  However,  it  has  been 
shown  that  the  photodetector  arrays  in  these  compact  imaging  systems  can  have  un¬ 
desirable  electronic  effects  that  degrade  range  precision.  A  key  undesirable  effect  that 
has  been  observed  in  photodetector  arrays  is  the  phenomenon  of  gain  variation. 

1.1.4  Airborne  LADAR  Systems.  Airborne  LADAR  systems  are  used  for 
three-dimensional  mapping  and  applications  such  as  target  recognition.  Typical  sys¬ 
tems  transmit  beams  through  apertures  that  are  several  inches  wide  and  use  pulses 
that  are  a  few  nanoseconds  in  duration.  If  the  pulse  reflects  off  a  target  that  is  tilted 
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Figure  1.2:  Block  diagram  describing  a  LADAR  system’s  re¬ 
ceiver  and  output  signal. 


with  respect  to  the  line-of-sight,  then  the  reflection  process  will  elongate  the  received 
signal  as  compared  to  transmitted  pulse.  If  the  range  is  more  than  a  few  kilometers 
and  the  target  is  tilted  more  than  about  forty-five  degrees,  the  increase  in  the  width  of 
the  received  pulse  produces  a  significant  drop  in  range  precision.  The  loss  in  precision 
due  to  target  tilt  is  quantified  in  this  dissertation. 

1.2  LADAR  Signals 

1.2.1  LADAR  Receivers.  A  LADAR  system  makes  measurements  of  target 
parameters  by  processing  the  laser  signal  detected  by  its  receiver.  A  LADAR  receiver 
block  diagram  is  shown  in  Fig.  1.2.  Transmitted  laser  light  reflects  off  of  objects  in  the 
scene.  Some  of  the  reflected  laser  light  is  intercepted  by  receiver  optics.  Those  optics 
focus  the  light  onto  a  photodetector,  which  is  small  compared  to  the  receiver  aperture. 
The  photodetector  converts  the  focused  optical  held  into  an  electrical  signal.  Typical 
detectors  are  photomultiplier  tubes  (PMT)s,  p-i-n  diodes,  and  avalanche  photodiodes 
(APD)s. 

There  is  noise  in  the  received  signal.  One  source  of  noise  is  background  light 
detected  by  the  photodetector.  The  background  light  level  is  roughly  constant  over 
the  time  it  takes  to  record  a  laser  pulse.  Thus,  it  interferes  with  detection  of  that 
pulse  and  degrades  measurements  of  signal  parameters  such  as  time  of  flight.  There 


5 


are  also  electronic  noise  sources  in  the  system.  Photodetectors  generate  false  signals 
called  dark  current.  The  signal  quantization  process  (analog-to-digital  conversion) 
also  introduces  error.  These  processes  are  all  sources  of  noise  in  the  received  signal  [30]. 

A  LADAR  must  make  multiple  target  range  measurements  to  create  a  three- 
dimensional  image.  Arrays  of  avalanche  photodiodes  are  used  in  LADAR  receivers  to 
create  three-dimensional  images.  Avalanche  photodiodes  convert  laser  light  into  an 
electrical  signal  called  the  photocurrent.  The  photocurrent  is  related  to  the  irradiance 
of  the  received  light  by  a  parameter  called  the  gain.  Ideally,  the  gain  would  remain 
constant  over  time.  However,  it  has  been  demonstrated  experimentally  that  the  gain 
varies  in  some  situations  due  to  limitations  in  photodiode  array  design.  It  has  also 
been  shown  that  gain  equalization  techniques  can  reduce  some  of  the  gain  variation 
in  the  data. 

In  order  to  characterize  gain  variation  in  LADAR  detector  arrays,  nonlinear 
models  are  also  being  used  to  describe  received  LADAR  signals.  The  nonlinear  models 
trade  simplicity  for  realism.  It  has  been  shown  in  a  real  LADAR  that  the  nonlinear 
effects  are  signihcant  and  that  attempts  should  be  made  to  correct  for  these  effects 
[37].  Correction  methods  were  proposed,  but  little  work  was  done  to  determine  range 
precision  in  the  nonlinear  signal  model. 

1.2.2  LADAR  Signals.  The  time  between  photon  arrivals  at  the  detector 
varies.  The  fluctuating  number  of  photons  arriving  in  a  time  interval  is  called  shot 
noise.  Shot  noise  is  a  source  of  error  in  LADAR  measurements.  Goodman  has  shown 
that  when  measuring  coherent  or  partially  coherent  light,  the  number  of  signal  pho¬ 
tocounts  observed  in  a  time  interval  is  a  negative  binomial  random  variable  [11, 13]. 
Negative  binomial  random  variables  have  Fano  factors  (the  ratio  of  the  variance  to  the 
mean  [7])  greater  than  one.  A  high  Fano  factor  indicates  that  the  signal  fluctuation 
is  large  compared  to  the  signal  energy.  LADAR  systems  are  often  designed  to  reduce 
the  Fano  factor  through  the  use  of  diversity.  Spatial  averaging  of  statistically  inde¬ 
pendent  laser  speckle  cells  and  the  detection  of  multiple  polarizations  are  common 
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methods  used  to  increase  diversity.  As  diversity  increases,  the  Fano  factor  approaches 
one,  and  the  negative  binomial  distribution  converges  to  the  Poisson  distribution.  Re¬ 
gardless  of  the  diversity,  the  Poisson  distribution  is  also  an  accurate  approximation 
to  the  negative  binomial  distribution  whenever  the  mean  number  of  photocounts  is 
low  [13].  Because  the  Poisson  distribution  can  accurately  approximate  the  negative 
binomial  distribution,  it  is  used  in  this  dissertation  to  model  the  number  of  observed 
photocounts  in  a  LADAR  signal. 

In  addition  to  the  received  laser  pulse,  signals  are  produced  by  noise  sources. 
Noise  sources  include  background  photons,  dark  current,  and  thermal  noise  [21]. 
These  sources  of  noise  corrupt  the  LADAR  signal  and  degrade  the  precision  of  range 
measurements.  Arrivals  of  background  photons,  dark  counts,  and  thermal  noise  counts 
are  random.  In  [21]  and  [31],  it  is  assumed  that  these  noise  sources  are  Poisson  pro¬ 
cesses.  This  model  is  reasonable  since  it  only  requires  independent  exponentially- 
distributed  times  between  noise  incidents.  The  number  of  counts  produced  by  a 
Poisson  process  within  a  specified  time  interval  is  a  Poisson  random  variable.  Sums 
of  statistically  independent  Poisson  random  variables  are  Poisson  with  mean  equal 
to  the  sum  of  the  means  of  the  constituent  random  variables.  Thus,  the  entire  ob¬ 
served  LADAR  signal  (laser  pulse  and  noise)  can  be  modeled  as  a  sequence  of  Poisson 
random  variables. 

The  Poisson  distribution  accurately  characterizes  the  signals  recorded  by  shot- 
noise  limited  LADAR  systems.  However,  the  noise  in  a  LADAR  signal  is  sometimes 
dominated  by  electronic  noise  sources  [9].  In  that  case,  the  stochastic  component  of 
the  LADAR  signal  can  be  modeled  with  a  Gaussian  distribution  [22].  The  Gaussian 
distribution  is  often  easier  to  work  with  than  the  Poisson  distribution. 

1.2.3  Signal  Parameter  Estimation.  LADAR  signal  processing  involves  es¬ 
timating  parameters  such  as  range,  amplitude,  and  background  level.  Methods  of  pa¬ 
rameter  estimation  include  maximum  likelihood  estimation  (MLE)  and  least  squares. 
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In  some  cases,  the  parameter  estimates  are  derived  empirically  because  there  are  often 
simple  and  nearly-optimal  estimators  of  signal  parameters  available. 

Though  many  different  parameters  can  be  estimated,  the  key  parameter  is  usu¬ 
ally  range  to  target.  In  order  to  make  a  precise  measurement  of  target  range,  the 
received  laser  pulse  is  usually  match-hltered.  After  hltering,  a  peak-htting  routine  is 
used  to  calculate  a  range  estimate  to  the  target.  However,  randomness  in  the  signal 
creates  fluctuations  in  the  output  of  the  peak-htting  algorithm.  These  huctuations 
cause  error  in  the  range  estimates.  The  distribution  of  the  observed  signal  is  as¬ 
sumed  known.  However,  the  combination  of  matched-hltering  and  peak-htting  makes 
it  difficult,  if  not  impossible,  to  derive  the  distribution  of  the  actual  range  estimate. 
Instead,  the  range  precision  is  analyzed  using  the  Cramer- Rao  Lower  Bound  (CRLB). 
The  CRLB  provides  a  lower  limit  on  the  variance  of  any  unbiased  estimate  of  a  pa¬ 
rameter  such  as  range.  The  CRLB  is  valid  regardless  of  the  exact  hltering  and  htting 
algorithms  used,  as  long  as  those  algorithms  produce  unbiased  range  estimates. 


II.  Problem  Background 

Typical  imaging  LADAR  systems  transmit  short  pulses  of  laser  light.  Modern  systems 
use  pulses  that  are  only  a  few  nanoseconds  in  duration.  The  pulses  scatter  off  of 
objects  in  a  scene.  The  backscattered  laser  light  is  detected  by  the  LADAR’s  receiver. 
The  received  signal  is  digitized  and  converted  to  an  estimate  of  range  to  the  target 
based  on  the  laser  pulse’s  time  of  flight  (TOF).  Multiple  estimates  of  target  range  can 
be  combined  into  a  three-dimensional  image  of  a  scene. 

Range  precision  is  limited  by  the  length  of  the  transmitted  laser  pulse,  the 
pulse’s  shape,  and  noise  sources  in  the  LADAR  system.  In  order  to  predict  limits 
on  LADAR  range  precision,  mathematical  models  of  the  shape  of  the  laser  pulse  are 
needed.  Statistical  methods  are  employed  to  simulate  the  noise  sources  present  in  the 
LADAR  system.  In  this  chapter,  laser  beam  models  used  in  the  literature  and  in  this 
dissertation  are  described.  The  stochastic  methods  used  to  simulate  noise  sources  are 
also  presented. 

Laser  beams  transmitted  by  LADAR  systems  are  characterized  by  their  shape 
transverse  to  the  line  of  sight  (LOS)  and  their  shape  along  the  line  of  sight.  The 
transverse  shape  describes  the  beam’s  irradiance  pattern.  As  the  beam  propagates, 
its  transverse  shape  varies  due  to  diffraction.  The  beam  is  also  described  by  a  shape 
along  the  LADAR  LOS.  This  shape  models  the  instantaneous  power  in  the  beam  and 
can  be  thought  of  as  a  temporal  function.  The  beam’s  shape  along  the  LOS  varies 
little  as  the  beam  propagates.  However,  the  pulse  shape  received  by  the  system  can 
be  very  different  from  the  transmitted  pulse  shape  because  of  the  reflection  process. 

Throughout  this  dissertation,  it  is  assumed  that  the  transverse  beam  shape  is 
Gaussian.  The  reasons  for  this  assumption  and  its  consequences  are  discussed  in 
Sec.  2.1.  Several  different  transmitted  temporal  (along  line  of  sight)  pulse  shape 
models  are  used  in  the  literature.  In  this  dissertation,  either  a  truncated  parabola 
or  a  Gaussian  function  is  used  to  describe  the  temporal  shape  of  a  laser  pulse.  The 
parabola  accurately  models  the  steep  leading  edge  of  real  laser  pulses.  It  also  proves 
easier  to  work  with  than  other  temporal  pulse  shape  models  when  deriving  bounds 
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on  range  precision  in  the  presence  of  Poisson-distribnted  noise.  The  parabolic  pulse 
model  is  discussed  in  Sec.  2.1.6.  The  Gaussian  shape  does  a  good  job  of  modeling  the 
smooth  trailing  edge  of  a  laser  pulse.  The  Gaussian  model  is  described  in  Ghapter 

IV. 

2.1  Laser  Beam  Propagation 

The  irradiance  pattern  of  the  laser  beam  transmitted  by  a  LADAR  system  can 
be  many  different  shapes.  The  most  common  shape  is  a  Gaussian  pattern.  However, 
top-hat  beams,  which  have  a  flat  irradiance  pattern  as  they  exit  the  LADAR  trans¬ 
mitter  aperture,  can  also  be  used.  There  are  also  higher-order  modes  described  by 
Hermite-Gaussian  functions  in  rectangular  coordinates  [42]  and  by  Laguerre- Gaussian 
functions  in  cylindrical  coordinates  [3].  The  fundamental  mode  created  by  a  laser 
resonator  is  the  called  the  zero-order  transverse  electromagnetic  or  TEMOO  mode. 
Typical  systems  transmit  TEMOO  Gaussian  beams. 

Zero-order  transverse  electromagnetic  Gaussian  beams  are  usually  used  in  LADAR 
systems  because  they  produce  the  best-quality  diffraction  patterns  in  the  far  field. 
Beam  quality  is  described  in  detail  in  the  following  section. 

2.1.1  Beam  Quality.  In  [39-41],  Siegman  developed  a  beam  quality  metric. 
The  metric  is  denoted  M^.  This  metric  is  the  product  of  the  beam’s  normalized 
irradiance  standard  deviation  and  the  beam’s  normalized  spatial  frequency  standard 
deviation  along  an  axis  transverse  to  the  LADAR  line  of  sight.  A  laser  beam  actually 
has  two  beam  quality  measurements.  There  is  one  measurement  along  each  axis 
perpendicular  to  the  LOS  (the  a:-axis  and  the  y-axis).  The  following  equations  used  to 
calculate  beam  quality  show  it  being  computed  along  the  x-axis.  It  is  straightforward 
to  modify  them  to  calculate  beam  quality  along  the  y-axis  in  case  that  measurement 
{My)  is  desired. 
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Beam  quality  along  the  x-axis  is  [39] 


Ml  ^ 


(2.1) 


The  normalized  spatial  variance  is  a  function  of  the  beam’s  irradiance  pattern 
Io{x,y).  Irradiance  has  units  of  power  per  area.  Irradiance  is  the  magnitude  squared 
of  optical  field  Uo{x,y). 

Io{x,y)  =  \Uo{x,y)f  .  (2.2) 


The  terms  irradiance  and  intensity  are  often  used  interchangeably  in  the  optics 
literature  [3].  However,  in  the  field  of  radiometry,  intensity  has  units  of  power  per 
area  per  solid  angle  [56].  Therefore,  intensity  could  also  be  defined  as  irradiance  per 
solid  angle.  To  avoid  confusion,  only  the  term  irradiance  is  used  in  this  dissertation 
when  referring  to  the  optical  held  Uo{x,y). 

The  normalized  spatial  variance  is 


OO  OO 

/  /  {x  -  xYlo{x,y)dxdy 

2  — OO  — OO 

OO  OO 

I  I  Ioix,y)dxdy 

—  OO  — OO 


(2.3) 


The  variable  x  in  Eq.  2.3  is  the  beam’s  irradiance  centroid  in  the  x-direction.  The 
beam  centroid,  which  is  the  beam’s  center  of  gravity,  is  (see  p.  662  of  [38]) 


OO  OO 

/  /  xIo{x,y)dxdy 

—  OO  — OO 


/  /  Io{x,y)dxdy 

—  OO  — OO 


(2.4) 


Calculation  of  the  normalized  spatial  frequency  variance  begins  with  com¬ 
putation  of  the  two-dimensional  Fourier  transform  of  the  optical  held  Uo{x,y).  This 
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function  is  denoted  UQ{sx,Sy).  This  Fourier  transform  is 


Uoi^Sxy  Sy^ 


Uo{x,  y)  exp  [— i27r  {sxX  +  Syy)]  dxdy 


(2.6) 


—  oo  — oo 


where  i  =  a/^.  The  spatial  frequency  distribution  Iq^Sx-,  Sy)  is  the  magnitude  squared 
of  this  Fourier  transform. 


d^oi^Sx^  Sy'j  Uoi^SxySy^ 


(2.6) 


The  normalized  variance  of  the  spatial  frequency  distribution  is 

OO  OO 

J"  J  Sy^dSxdSy 

2  — oo  — oo 

^Sx  oo  oo  ^ 

J"  J  Sy^dSxdSy 

—  oo  — oo 


(2.7) 


where  Sx  is  the  centroid  of  the  spatial  frequency  distribution.  This  variable  is  com¬ 
puted  in  the  same  way  the  irradiance  centroid  was  calculated  in  Eq.  2.4. 


OO  OO 

J  J  ^xloi^^X-)  Sy^dSxdSy 

_  — OO  — oo 

S'T  - 

^  oo  CO 

J"  J  Sy^dSxdSy 

—  OO  — oo 


(2.8) 


Siegman  showed  in  [39]  that  for  any  arbitrary  laser  beam 


Mx  >  1 


(2.9) 


with  equality  if  and  only  if  the  irradiance  pattern  has  a  Gaussian  shape.  Thus,  the 
beams  with  the  best  quality  have  Gaussian  irradiance  patterns  and  Gaussian  beams 
are  almost  always  used  in  LADAR  systems. 
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Combining  Eqs.  2.1  and  2.9  shows  that 


(2.10) 

This  bound  on  the  product  of  space  and  spatial  frequency  bandwidths  is  equivalent 
to  the  bound  on  time  and  frequency  that  appears  in  literature  on  radar  and  signal 
processing  [51]. 


2.LS  Fresnel  Propagation  of  Gaussian  Beams.  The  optical  held  after  prop¬ 
agating  distance  z  is  denoted  [/o(r,  z).  The  held  for  a  Gaussian  beam  at  the  aperture 
{z  =  0)  is  given  by  Eq.  8  of  [3]. 


Go(r,0) 


(2.11) 


where 


a  = 


kW^ 


i 


(2.12) 


In  Eq.  2.12,  Wq  is  the  beam  size,  Fq  is  the  focus  distance,  and  k=  2tt/\  is  the  wave 
number.  The  beam  size  (ITo)  is  the  e~^  radius  of  the  magnitude  of  the  held  at  the 
aperture.  In  other  words,  it  is  the  e~^  radius  of  |Go(r,  0)|.  The  power  in  this  beam  is 
denoted  P. 

27r  CO  OO 

P  =  j  j \Uoir,0)\^rdrde  =  2%  j  |Go(r,0)|Vdr  (2.13) 

0  0  0 

where  the  radial  symmetry  in  the  held  has  been  exploited.  The  power  in  Eq.  2.13 
is  measured  at  an  instant.  In  LADAR  systems  transmitting  laser  pulses,  the  power 
varies  with  time. 

The  Fresnel  propagation  of  the  held  from  Eq.  2.11  is 


Go(r,z)  = 


1/2 


2P 


t^Wq  j  1  -|-  iaz 


exp 


ikz 


akP 


2(1  -|-  iaz) 


(2.14) 


13 


This  is  Eq.  24  of  [3].  The  irradiance  after  propagating  distance  z  is 


h{r,z) 


2P/'kW^ 


(1  -  z/FoY  +  (Az/ttIT^)- 


exp 


2(r/fTo) 


(1  -  z/FoY  +  {Xzl-nWiY 


(2.15) 


When  the  beam  is  collimated  {Fq  — >  cx)),  the  propagated  irradiance  is 


2F/ttW^ 


2^2 


1  +  (Az/vrlTo^) 


exp 


2  (r/Wo 


2^2 


1  +  {Xz/nW^) 


(2.16) 


When  the  beam  is  at  the  focus  distance  {z  =  Fq),  the  propagated  irradiance  is 


Io{r,Fo)  = 


2F  /nW, 


ttW^ 


2\  2 


AFn 


exp 


V  AFo  j 


27iFW^ 


exp 


J^Wor\ 

\XF^  ) 


(2.17) 


2.1.3  Gaussian  Beam  Size.  The  beam  size  ITo  is  measured  at  the  LADAR 
system’s  transmitter  aperture.  As  the  beam  propagates,  its  size  varies.  The  size  as  a 
function  of  range  is  calculated  from  the  held  magnitude  equation,  which  is  the  square 
root  of  the  beam’s  irradiance.  The  Fresnel  propagated  held  magnitude  is 


mr,z)\  (2.18) 


2F/7iW^ 

1/2 

{r/Wof 

[{l-z/FY^  +  {Xz/T:Wif\ 

{l-z/Fof  +  {\z/'KW^f_ 

The  beam  size  at  range  z,  which  is  denoted  IW,  is 


W,  =  Wo 


(2.19) 


This  is  the  radius  at  which  the  held  magnitude  is  below  the  peak  magnitude  by  a 
factor  of  e~^. 
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The  beam  size  at  range  2;  when  the  beam  has  been  focused  is 


AFq 

vrlTn 


(2.20) 


When  the  beam  is  collimated  {Fq  — >■  cxo),  the  beam  size  at  range  2;  is  (see  also  Eq. 
5.85  of  [28]) 


1  + 


Az 


1/2 


(2.21) 


The  equation  for  the  propagated  beam  size  Wz  (Eq.  2.19)  can  be  used  to  write 
the  field  magnitude  and  irradiance  equations  in  simpler  forms.  The  field  magnitude 
from  Eq.  2.18  can  be  written 


\Uo{r,z)\  = 


2P 


7rlT'2 


exp 


(2.22) 


The  propagated  irradiance  from  Eq.  2.17  can  be  expressed 


/  ,  2P  /  r2 


(2.23) 


2.1.4  Beam  Divergence  Angle.  In  the  far  held,  the  beam’s  size  is  approxi¬ 
mately  proportional  to  the  distance  of  propagation.  Because  of  this  linear  dependence, 
the  beam’s  size  is  often  measured  in  terms  of  the  angle  between  the  axis  of  propaga¬ 
tion  and  the  e~^  held  magnitude.  The  TEMOO  Gaussian  beam  divergence  angle  6^  is 


(Eq.  5.86  of  [28]) 


0ci  = 


A 


(2.24) 


The  Gaussian  beam  divergence  angle  can  be  compared  to  the  divergence  an¬ 
gle  for  a  plane  wave  truncated  by  an  aperture  with  diameter  Da.  The  plane  wave 
divergence  angle  9p  is  (Eq.  5.87  of  [28]) 


e„  =  1.22  — . 

^  D 


(2.25) 
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Field  Magnitude 


Field  Magnitude 


Figure  2.1:  (a)  Propagation  of  a  4  cm  beam  that’s  focused  at  2  km. 

(b)  Propagation  of  a  4  cm  beam  that  is  focused  at  inhnity  (collimated). 

The  scaling  of  1.22  in  the  previous  equation  is  the  first  zero  of  Ji(7ra:)  where  Ji  is  the 
modihed  Bessel  function  of  the  first  kind  of  order  one.  (See  Table  4.1  of  [13].) 

2.1.5  Propagation  Examples.  Examples  of  propagated  TEMOO  Gaussian 
beams  are  shown  in  Fig.  2.1.  In  Fig.  2.1a,  the  beam  is  focused  at  2  km.  The  beam  size 
is  4  cm  and  the  wavelength  is  1.6/im.  Notice  that  the  beam  is  narrower  after  1  km  of 
propagation  than  at  the  2  km  focus  distance.  A  beam  focused  at  a  particular  distance 
has  the  narrowest  possible  value  at  that  distance.  But  it  could  be  even  narrower  at  a 
shorter  distance  than  the  focus.  In  Fig.  2.1b,  propagated  collimated  Gaussian  beams 
are  shown.  In  the  collimated  case,  the  beam’s  size  increases  monotonically  with  range. 

The  Gaussian  beam  size  as  a  function  of  range  is  shown  in  Fig.  2.2a.  Beams 
focused  at  2,  5,  and  10  km  are  shown,  along  with  a  collimated  beam.  The  divergence 
angle  from  Eq.  2.24  is  shown  for  comparison.  Notice  that  the  beam  focused  at  short 
range  (2  km)  begins  to  diverge  faster  than  the  collimated  beam  at  longer  ranges.  If 
the  beams  propagate  far  enough,  any  Gaussian  beam  not  collimated  will  eventually 
begin  to  diverge  faster  than  the  Eq.  2.24  divergence  angle. 
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Gaussian  Beam  Divergence 


(a) 


Plane  Wave  Divergence 


(b) 


Figure  2.2:  (a)  Beam  divergence  for  Gaussian  beams, 

(b)  Beam  divergence  for  plane  wave  beams. 


The  plane  wave  divergence  angle  from  Eq.  2.25  is  shown  in  Fig.  2.2b.  The 
angle  for  10,  15,  and  20  cm  apertures  is  illustrated.  Note  that  the  15  cm  plane 
wave  divergence  angle  is  nearly  equal  to  the  Gaussian  beam  divergence  angle  for  a 
transmitted  beam  size  of  4  cm. 


2.1.6  Temporal  Laser  Pulse  Shape.  Derivation  of  the  Gramer-Rao  lower 
bound  on  range  precision  requires  a  mathematical  model  of  the  received  laser  pulse’s 
temporal  shape.  Several  different  LADAR  pulse  shape  models  are  used  in  the  litera¬ 
ture.  Both  a  Gaussian  shape  [22]  and  heavy-tailed  (asymmetric)  [15,45]  models  have 
been  previously  employed.  In  this  dissertation,  parabolas  and  Gaussian  shapes  are 
used  to  describe  the  temporal  shape  of  the  transmitted  laser  pulse. 

The  Gaussian  model  does  not  work  well  when  it  is  assumed  that  the  noise  is 
Poisson.  In  that  case,  GRLBs  are  very  difficult  to  calculate.  A  very  similar  problem 
was  considered  by  Winick  [55].  That  author  considered  the  problem  of  calculating 
the  GRLB  on  position  of  a  Poisson-distributed  Gaussian  shape  in  two- dimensions.  A 
GRLB  was  presented,  but  the  author  did  not  show  a  closed-form  solution.  However, 
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it  is  straightforward  to  calculate  the  range  CRLB  for  a  Gaussian  pulse  using  the  same 
sort  of  numeric  expressions  that  are  in  [4, 55] . 

A  solution  for  the  shape  of  a  pulse  created  by  a  Q-switched  laser  is  presented 
by  Siegman  in  [38].  However,  that  solution  is  expressed  as  an  integral  that  must  be 
evaluated  numerically.  Cramer-Rao  lower  bound  analysis  on  such  a  signal  would  be 
difficult.  Even  if  an  exact  solution  for  the  shape  of  the  transmitted  laser  pulse  was 
used,  the  reflection  process  can  distort  the  shape  of  that  pulse,  causing  the  received 
pulse  to  have  a  different  shape.  The  truncated  parabola  used  in  this  dissertation  is 
not  meant  to  be  an  exact  representation  of  the  true  received  signal.  Rather,  it  is  an 
adequate  approximation  of  received  pulse  shape  and  it  produces  a  relatively  simple 
range  CRLB  in  the  Poisson  noise  case. 

The  various  laser  pulse  shape  models  each  have  advantages  and  disadvantages. 
The  asymmetric  models  from  [15, 45]  accurately  simulate  the  relatively  fast  rise  of  a 
laser  pulse.  However,  those  models  can  overestimate  the  amount  of  pulse  energy  in 
the  tail.  The  asymmetry  also  makes  analytic  analysis  of  range  precision  more  difficult. 

The  parabolic  shape  generally  models  the  fast  rise  of  a  real  laser  pulse  better 
than  the  Gaussian  model.  However,  the  Gaussian  model  approximates  the  tail  of  the 
laser  pulse  better  than  the  parabolic  model.  The  parabolic  model  is  truncated  and 
does  not  have  the  tail  that  is  present  in  a  laser  pulse.  But  since  only  a  small  fraction 
of  the  pulse  energy  is  present  in  the  tail  of  a  real  laser  pulse,  the  inverted  parabola  is 
a  reasonable  approximation. 

2.2  LADAR  Detectors  and  Signals 

2.2.1  Photodetectors  and  their  Output  Signals.  A  portion  of  the  propagated 
field  is  reflected  off  the  target.  Some  of  the  reflected  field  is  intercepted  by  the  LADAR 
system’s  receiver  optics,  which  focus  that  field  onto  a  photodetector.  A  LADAR’s 
detector  converts  the  received  optical  power  to  a  photocurrent. 
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The  photocurrent  is  denoted  i{t).  It  is  a  function  of  the  optical  held  and  the 
detector’s  responsivity. 

OO  OO 

^(^)  =  /  J  Ad{x,y)^{x,y)\Uo{x,y)\‘^dxdy  (2.26) 

—  CO  — OO 

where  Ad{x,y),  written  as  a  function  of  spatial  coordinates  x  and  y,  is  the  binary 
detector  area  function.  This  function  represents  the  detector’s  active  region.  In  this 
case,  the  optical  held  UQ{x,y)  is  evaluated  in  the  detector  plane. 

The  variable  in  Eq.  2.26  is  detector  responsivity.  The  responsivity  is 

given  by  [13] 

»(..,)  =  (2.27) 

where  G  is  the  detector  gain,  q  is  the  electron  charge  (1.602  x  10“^®  Coulombs),  r]{x,  y) 
is  the  photodetector’s  quantum  efficiency,  h  is  Planck’s  constant  (6.626  x  10“^^  Joule- 
seconds),  and  V  is  the  optical  frequency  of  the  laser.  The  units  of  responsivity  are 
electric  charge  per  unit  of  energy  (e.g..  Coulombs  per  Joule),  or  equivalently,  current 
per  power  (e.g..  Amps  per  Watt).  The  spatial  coordinates  are  frequently  dropped 
from  the  responsivity  and  quantum  efficiency  terms  because  these  values  are  often 
nearly  uniform  over  the  photodetector. 

Targets  observed  by  imaging  systems  are  usually  extended.  That  is,  they  are 
large  compared  to  the  size  of  the  laser  beam.  When  extended  targets  are  observed, 
none  of  the  transmitted  laser  power  misses  the  target.  However,  some  laser  power 
is  absorbed  by  the  target  materials.  The  materials  illuminated  by  LADAR  systems 
are  often  assumed  to  be  Lambertian.  This  means  that  reflected  laser  light  scatters 
equally  in  all  directions. 

The  power  received  Pr{t)  by  a  LADAR  system  illuminating  an  extended  Lam¬ 
bertian  target  can  be  approximated  as  [21,31] 

PAt)  =  '-^P,  (t  -  IRjc)  (2.28) 
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where  Pt{t)  is  the  power  transmitted  at  time  t,  r]s  is  the  system  efficiency,  p  is  the 
target’s  reflectivity,  is  the  area  of  the  receiver,  R  is  the  range  to  the  target,  and 
c  is  the  speed  of  light.  The  system  efficiency  rjs  is  the  product  of  the  various  loss 
factors  such  as  atmospheric  transmission  and  optical  transmission.  In  a  system  using 
detector  arrays,  the  received  power  is  divided  between  the  pixels  in  the  array. 

2.2.2  LADAR  Signal  Model.  The  following  model  was  used  in  Flash  light 
detection  and  ranging  precision  limits  for  returns  from  single  opague  surfaces  via 
Cramer-Rao  bounds”  by  Cain,  et  ah  [4]  to  represent  the  mean  number  of  photoelec¬ 
trons  per  sample  I{tk,  Xn,ym)  in  an  imaging  LADAR  system.  The  model  describes 
the  received  laser  pulse  as  an  upside-down  parabola.  This  LADAR  signal  model  is 

I{tk,Xn,ym)  (2.29) 

=  G(x„.  y„)  1 1  -  I  rect 

-\-B  i^Xnj  ym)  T  Qsitk'!  Xni  ym) 

where  gain  term  G{xn,ym)  is  the  peak  mean  number  of  photoelectrons  per  sample 
from  the  target,  Rrt{xn,ym)  is  the  round-trip  range  to  the  target,  t^  is  the  time  at 
which  the  k-th  sample  is  collected,  is  the  pulse  width,  c  is  the  speed  of  light,  and 
B(xn,ym)  represents  a  bias  level  on  each  pixel  in  the  LADAR’s  detector  array. 

The  variables  m  and  n  in  Eq.  2.29  are  the  spatial  indices  in  the  LADAR  image. 
These  indices  specify  the  coordinate  of  the  detector  in  the  receiver’s  detector  array. 
They  index  x^  and  |/„,  which  are  the  spatial  coordinates.  Since  much  analysis  in 
this  dissertation  (e.g.,  CRLB  derivations)  is  done  on  only  the  data  from  one  pixel  at 
a  time,  the  spatial  coordinates  indexed  by  these  variables  are  sometimes  dropped  in 
equations.  However,  it  has  been  noted  that  a  large  signal  incident  on  one  detector 
can  effect  adjacent  detectors  in  an  imaging  LADAR’s  focal  plane  array  [37].  In  such 
a  case,  it  is  not  valid  to  ignore  data  from  the  other  detectors  in  the  array,  and  the 
spatial  coordinates  must  be  included  in  any  analysis. 


Rrt{Xniym)  tfcC 
2cp^ 
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The  bias  term  is  the  mean  number  of  noise  photoelectrons  per  sample  observed 
by  the  detectors  in  the  LADAR  system.  It  represents  the  mean  number  of  noise 
counts  due  to  noise  sources  such  as  background  light  and  dark  current.  The  term 
qs(tk,Xn,ym)  is  a  stochastic  noise  component  in  each  voxel  of  the  three-dimensional 
image.  Inclusion  of  this  term  allows  the  noise  level  to  vary  for  different  samples  from 
a  particular  pixel.  The  units  of  qsitk,Xn,ym)  are  photoelectrons  per  sample.  The 
stochastic  noise  term  qs{tk,  Xn,  ym)  is  sometimes  dropped  from  this  equation  since  the 
gain  and  bias  terms  may  accurately  model  the  signal  by  themselves. 

The  rect(x)  function  in  Eq.  2.29  is  the  standard  rectangle  function  from  the 
optics  and  signal  processing  literature.  It  is  defined  in  [13]. 


( 

1 

when 

|a;|  <  1/2 

rect(a:)  =  < 

1 

2 

when 

|a;|  =  1/2  . 

(2.30) 

0 

when 

|a;|  >  1/2 

The  parabolic  pulse  shape  used  in  the  signal  model  from  Eq.  2.29  is  illustrated 
in  Fig.  2.3.  This  illustration  shows  that  the  pulse  width  p^,  in  the  parabolic  signal 
model  is  the  distance  between  the  points  that  are  three-quarters  of  the  maximum  pulse 
height.  The  entire  pulse  width  is  actually  2p^.  The  bias  illustrated  in  Fig.  2.3  is  one 
{B  =  1)  and  the  gain  is  three  (G  =  3).  The  round-trip  range  is  ten  meters  {Rrt  =  10 
m).  It  has  been  assumed  that  the  stochastic  noise  term  is  zero,  i.e.,  qs(tk,Xn,ym)  =  0, 
in  this  hgure. 

2.2.3  Photoelectron  Probability  Density  Functions.  The  photocurrent  is  pro¬ 
portional  to  the  number  of  photoelectrons  observed  by  the  LADAR  photodetector. 
The  observed  number  of  photoelectrons  is  denoted  D{tk,  Xn,  ym)-  The  process  of  con¬ 
verting  a  photon  to  a  photoelectron  is  called  a  photoevent.  The  number  of  photoevents 
that  occur  in  some  time  duration  is  called  the  number  of  photocounts.  These  terms 
photoevent  and  photocounts  are  described  in  detail  in  Sec.  9.1  of  Statistical  Optics  by 
Goodman  [12]. 
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Parabolic  Pulse  LADAR  Signal  Model 


Figure  2.3:  Illustration  of  the  parabolic  pulse  LADAR  signal 
model.  This  model  is  used  for  estimator  derivations  and  for 
CRLB  calculations. 
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Through  much  of  this  dissertation,  it  is  assumed  that  the  number  of  photocounts 
D{tk,Xn,ym)  has  a  Poisson  distribution  with  mean  I{tk,Xn,ym)  given  by  Eq.  2.29. 


It  has  been  shown  that  the  integrated  irradiance  of  coherent  light  reflected  off  of 


a  diffuse  surface  actually  has  a  negative  binomial  distribution  (see  Sec.  A. 2).  The 
number  of  photoelectrons  observed  by  a  LADAR  system’s  receiver  D{tk,Xn,ym)  is 
proportional  to  this  integrated  irradiance  (see  Sec.  9.2.3  of  [12]).  They  are  related  by 
the  gain  and  the  quantum  efficiency  of  the  detector  (often  a  photodiode  or  PMT)  in 
the  LADAR’s  receiver.  The  negative  binomial  distribution  is  parameterized  by  the 
mean  number  of  photoelectrons  and  the  diversity  parameter  M..  If  the  mean  number 
of  photoelectrons  is  fixed  and  the  diversity  M.  becomes  large  in  the  negative  binomial 
distribution,  then  that  distribution  converges  to  the  Poisson  distribution  [12].  It  has 
also  been  shown  that  as  the  mean  number  of  photoelectrons  approaches  zero,  that 
the  negative  binomial  distribution  is  approximately  equal  to  the  Poisson  distribution 
(see  Sec.  A. 3). 

Typical  direct  detection  LADAR  systems  have  a  relatively  high  diversity.  Thus, 
the  Poisson  distribution  is  often  used  as  an  accurate  approximation  to  the  true  nega¬ 
tive  binomial  distribution  that  describes  the  number  of  photoelectrons  in  the  LADAR 
signal.  For  example,  the  system  discussed  in  [4]  had  a  diversity  of  about  nine  hundred 
{M.  ~  900).  The  Poisson  distribution,  used  to  model  the  number  of  photoelectrons 
counted  in  some  time  interval,  is 


Prob  [D(4,  Xn,  ym)  =  d{tk,  Xn,  ym)]  = 


(2.31) 


where 

d(tky  ym)  0,  1,  2,  ...  . 


The  variable  d{tk,  Xn,  ym)  is  a  realization  of  photocurrent  random  variable  D(tk,  Xm  ym)- 
Variable  I{tk,Xn,ym)  is  the  mean  value  for  each  sample. 

2.2.4  Single  Detector  LADAR  Signal  Model.  At  this  point,  the  mean  pho¬ 
tocurrent  Eq.  2.29  is  simplihed.  First,  it  is  assumed  that  only  one  pixel  in  the  LADAR 
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detector  array  is  being  analyzed.  Thus,  the  spatial  indexing  coordinates  {xm,  Vn)  from 
Eq.  2.29  have  been  dropped  in  order  to  simplify  the  notation.  Second,  it  is  assumed 
that  the  stochastic  noise  term  qs(tk,Xn,ym)  is  zero. 

The  parabolic  pulse  LADAR  signal  model  without  spatial  dependence  is 


(2.32) 


This  equation  is  a  simplihcation  of  the  model  from  [4]  (which  is  also  shown  in  Eq.  2.29 


of  this  dissertation).  The  more  general  model  includes  spatial  dependence  (denoted 


by  variables  Xn  and  i/m)  and  a  stochastic  noise  term  qsitk,Xn,ym)  for  each  voxel  in 
the  three-dimensional  image. 

In  this  chapter,  it  is  assumed  that  each  signal  received  by  the  LADAR  system 
is  statistically  independent  from  all  others.  Therefore,  the  each  pulse  can  be  treated 
individually  for  purposes  for  parameter  estimation  and  CRLB  derivations.  The  fact 
that  each  pulse’s  signal  is  independent  justihes  the  dropping  of  the  spatial  variables 
when  going  from  Eq.  2.29  to  Eq.  2.32. 

The  number  of  temporal  samples  recorded  by  the  (m,n)-th  detector  is  denoted 
K.  If  it  is  further  assumed  that  the  samples  are  statistically  independent  of  each 
other,  then  the  joint  PDF  of  the  photoelectrons  observed  at  each  sample  is  simply 
the  product  of  the  K  PDFs  for  each  individual  sample.  The  assumption  of  statistical 
independence  between  samples  comes  from  the  fact  that  the  photocurrent  is  assumed 
to  be  a  Poisson  process.  Poisson  processes  are  memoryless.  Therefore,  any  sample 


24 


has  no  dependence  on  all  previons  samples.  The  joint  PDF  is 


Prob  [D{ti)  =  d{ti),  D{t2)  =  ^(^2),  •  •  • ,  Di^x)  =  d(tx)] 


=  p[d{ti),d{t2),. . .  ,d{tK)] 


=  n 


/(tfc)'^hfc)e  idk) 

d(4)! 


exp 

-  Htk) 

k=l 

n 

k=l 

K 

n  d{tky. 

k=l 


(2.33) 


where  d{tk)  is  the  observed  realization  of  the  Poisson  process  and  I{tk)  is  the  mean 
nnmber  of  photons  from  Eq.  2.32. 


2.3  LADAR  Range  Precision 

The  precision  of  range  estimates  in  LADAR  systems  has  been  stndied  in  previons 
pnblications.  This  precision  is  bonnded  nsing  Cramer-Rao  analysis.  Derivation  of  the 
CRLB  reqnires  a  statistical  model  for  the  received  signal.  That  model  inclndes  the 
noise  terms  present  in  the  received  signal.  Unknown  parameters  in  the  model,  snch 
as  amplitnde  of  the  received  signal  and  range  to  target,  must  be  estimated  from  the 
noisy  data.  The  model  is  used  in  the  Cramer-Rao  lower  bound  to  determine  precision 
limits  on  unbiased  estimates  of  the  unknown  parameters. 

2.3.1  Cramer-Rao  Lower  Bounds.  At  this  point,  the  Cramer-Rao  lower 
bounds  for  Eq.  2.32  are  calculated  using  the  distribution  from  Eq.  2.33.  In  the  K 
sample  LADAR  signal,  there  are  three  unknown  parameters  that  can  be  estimated. 
The  round-trip  range  to  target  Rrt,  the  gain  G,  and  the  bias  B  are  all  unknown 
parameters.  The  precision  of  estimates  of  the  unknown  parameters  can  be  analyzed 
using  the  Cramer-Rao  lower  bound.  The  CRLB  is  found  by  hrst  evaluating  the 
Fisher  information  matrix  (FIM).  The  Fisher  information  matrix  3{Rrt,G,B)  for 
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this  problem  is 


(2.34) 


where  l{Rrt,G,  B)  is  the  log-likelihood  function. 

The  likelihood  function,  which  is  denoted  L{Rrt,  G,B),  is  closely  related  to  the 
PDF  of  the  data.  It  is  the  PDF  written  as  a  function  of  the  unknown  parameters. 
The  PDF,  on  the  other  hand,  is  a  function  of  the  data  realization.  The  likelihood 
function  derived  from  Eq.  2.33  is  a  function  of  the  round-trip  range,  gain,  and  bias  in 
the  signal.  Specihcally,  the  likelihood  function  associated  with  the  distribution  from 
Eq.  2.31  is 
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(2.35) 


The  log-likelihood  function,  which  is  denoted  l{Rrt,G,  B),  is  simply  the  natural  log¬ 
arithm  of  the  likelihood  function.  The  log-likelihood  function  for  Eq.  2.31  is 


I  [Rrt,  G,  B\d{ti),  d{t2), . . . ,  d{tK)] 


(2.36) 
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The  second  derivatives  of  the  log-likelihood  function  with  respect  to  unknown  param¬ 
eters  Rrti  G,  and  B  are 


dH{Rrt,G,B) 
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(2.37) 
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The  elements  of  the  Fisher  information  matrix  are  the  negative  expectations  of  these 
second  derivatives 
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The  derivations  of  the  previous  equations  made  use  of  the  fact  that  the  mean  of  a 
Poisson  random  variable  Ditk)  is  litk).  That  is,  E[T)(4)]  =  litk).  Eqs.  2.43  -  2.48 
were  first  derived  in  [4]. 

The  functions  being  summed  in  Eqs.  2.43  -  2.48  are  illustrated  in  Fig.  2.4. 
The  round-trip  range,  gain,  and  bias  values  used  to  create  these  plots  are  the  same 
as  are  used  in  Fig.  2.3.  It  is  shown  in  Chapter  III  that  the  summations  in  these 
equations  can  be  approximated  using  integrals.  The  functions  in  Eqs.  2.44  and  2.45 
have  odd  symmetry  around  the  center  of  the  pulse.  Therefore,  those  equations  are 
approximately  zero. 

The  duration  of  the  pulse  is  2p^„.  The  sampling  duration,  which  is  denoted 

is 

td  =  K/f,  =  KAt  (2.49) 


where  fs  is  the  sampling  frequency.  The  sampling  frequency  is  the  reciprocal  of  the 
time  between  samples  in  the  system  (/„  =  1/At)  where  At  denotes  the  time  between 
samples. 


The  pulse  in  Eq.  2.32  is  centered  at  time  =  Rh/c.  Since  we  have  assumed  that 
the  entire  pulse  is  sampled,  it  is  required  that  <  Rrt/c  <  td—Pw  If  this  requirement 
is  not  satisfied,  then  the  pulse  is  not  completely  sampled  and  some  received  signal 
energy  is  missed. 

2.4  Gain  Variation  in  APD  Arrays 

In  Sec.  2.3,  a  model  was  presented  to  describe  the  received  signal  in  the  pho¬ 
todetector  array  of  an  imaging  LADAR.  The  model  from  Eq.  2.29  assumes  that  the 
gain  term  G  is  constant  across  the  entire  detector  array.  That  model  also  assumes 
that  the  gain  is  fixed  for  each  pixel  over  its  range  record.  Fixed  detector  gain  was 
also  assumed  in  [4]. 

However,  recent  research  has  demonstrated  that  there  can  be  gain  variations 
in  LADAR  detector  arrays  [37].  Gain  can  vary  for  a  particular  detector  within  an 
array  over  that  detector’s  sampling  time.  Gain  can  also  fluctuate  across  an  array  of 
photodetectors  at  any  particular  sample.  These  fluctuations  in  gain  were  observed  in 
the  Advanced  Scientific  Goncepts  LADAR  owned  by  AFRL  Sensors  Directorate.  The 
detector  used  in  that  system  is  an  array  of  avalanche  photodiodes. 

An  empirical  model  for  the  gain  in  an  avalanche  diode  was  developed  by  Miller 
[29].  Miller’s  model  is 


where  Vbias  is  the  bias  voltage,  Vm  is  the  breakdown  voltage,  and  p  is  an  empirical 
constant  used  to  model  the  device.  Some  imaging  LADARs  use  detector  arrays  that 
share  a  common  voltage  regulator,  which  is  the  component  that  sets  the  bias  voltage 
Vbias-  Sharing  the  voltage  regulator  saves  on  weight  and  power  in  the  system.  However, 
when  a  signal  is  received,  the  regulator’s  bias  voltage  drifts.  Fluctuations  in  the  bias 
voltage  cause  variations  in  the  gain  equation  (Eq.  2.50).  And  if  the  same  regulator  is 
used  by  all  detectors  in  the  system,  then  the  gain  will  vary  in  all  pixels. 
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Figure  2.4:  (a)  The  round-trip  range  second  derivative,  (b)  The  range-gain  deriva¬ 

tive.  (c)  The  range-bias  derivative,  (d)  The  gain  second  derivative,  (e)  The  gain-bias 
derivative,  (f)  The  bias  second  derivative. 
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Methods  for  compensation  of  gain  variation  were  proposed  and  tested  in  [37] ,  but 
further  work  was  also  recommended.  In  this  dissertation,  more  LADAR  data  collected 
with  an  ASC  system  is  analyzed.  The  gain  variation  is  observed  and  its  effect  on  range 
measurements  is  demonstrated.  Different  methods  of  gain  equalization  are  proposed 
and  tested  in  an  attempt  to  correct  for  the  gain  variation  in  the  avalanche  photodiode 
array. 

2. 5  Overview  of  Dissertation 

The  research  in  this  dissertation  begins  in  Chapter  III  with  a  Poisson-distributed 
signal  model  that  simulates  a  laser  pulse  received  by  a  shot-noise  limited  LADAR 
system.  This  signal  model  was  previously  analyzed  by  Cain  et  al.  in  [4].  In  this 
dissertation,  the  Cramer-Rao  lower  bound  analysis  work  from  that  paper  is  extended. 
Numerical  solutions  for  range,  gain,  and  background  level  CRLBs  were  formerly  pub¬ 
lished  in  that  paper.  In  this  dissertation,  new  closed-form  analytic  expressions  are 
derived  for  those  Cramer-Rao  lower  bounds. 

The  shot-noise  limited  signal  model  from  [4]  is  also  extended  to  cases  where 
there  are  multiple  received  signals.  These  cases  simulate  the  signal  in  a  LADAR 
observing  a  target  partially  obscured  by  forest  canopy.  Earlier  work  with  this  signal 
model  assumed  that  the  received  pulse  width  was  known.  However,  when  a  target  has 
depth  that  is  signihcant  compared  to  the  laser  pulse  width,  the  received  pulse  width 
is  larger  than  the  transmitted  width.  In  this  case,  the  known  width  assumption  is  not 
realistic.  In  Chapter  III,  the  effect  of  an  unknown  pulse  width  is  analyzed  and  the 
effect  on  signal  parameter  CRLBs  is  calculated. 

In  Chapter  IV  of  this  dissertation,  the  range  CRLB  for  a  Poisson-distributed 
signal  is  compared  to  the  CRLB  for  a  signal  mixed  with  Gaussian  noise.  The  compar¬ 
ison  assumes  the  background  Poisson  and  Gaussian  noise  variances  are  equal.  Proof 
that  the  Gaussian  noise  range  GRLB  is  a  lower  bound  for  the  Poisson  noise  GRLB  is 
presented.  Also,  the  effect  of  using  multiple  pulses  on  range  precision  of  a  shot-noise 
limited  LADAR  is  studied.  It  is  shown  that  for  a  laser  capable  of  dividing  a  hxed 
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amount  of  energy  into  multiple  pulses,  the  best  range  precision  is  achieved  when  all 
energy  is  transmitted  in  a  single  pulse.  Equivalently,  distributing  a  hxed  amount  of 
energy  into  multiple  pulses  increases  the  range  CRLB,  which  indicates  that  there  is 
a  decrease  in  precision. 

It  is  also  noted  in  Chapter  IV  that  there  are  cases  where  the  effect  of  target 
surface  orientation  on  range  precision  can  not  be  ignored.  These  cases  usually  involve 
systems  with  modest  aperture  sizes  (a  few  inches  in  diameter)  and  target  ranges  of  at 
least  several  kilometers.  However,  such  scenarios  are  common  for  airborne  LADAR 
systems.  If  the  target  surface  is  tilted  with  respect  to  the  LADAR  line  of  sight, 
then  there  is  a  drop  in  range  precision.  The  drop  in  precision  is  quantified  in  this 
dissertation. 

In  Chapter  V,  it  is  noted  that  detector  arrays  composed  of  avalanche  photodi¬ 
odes  (APDs)  are  being  used  in  current  imaging  LADAR  systems.  Current  arrays  have 
tens  of  thousands  of  APDs.  The  compact  size  and  high  numbers  of  pixels  makes  the 
use  of  APD  arrays  desirable  in  imaging  LADARs.  However,  gain  variation  in  APD 
arrays  can  produce  undesirable  effects.  Gain  variation  distorts  the  shape  of  received 
signals  and  introduces  bias  into  range  measurements. 

In  Chapter  V,  methods  of  gain  variation  correction  (gain  equalization)  are  pro¬ 
posed  and  tested.  It  is  shown  that  background  data  collected  by  an  APD  array 
can  be  used  to  decrease  the  gain  variation  and  improve  the  range  precision  in  three- 
dimensional  LADAR  images. 

2.6  Chapter  Summary 

Understanding  of  LADAR  systems  requires  knowledge  of  several  key  subjects. 
The  propagation  of  laser  light  was  described  in  Sec.  2.1  of  this  chapter.  Photodetectors 
are  used  to  convert  a  received  optical  held  to  an  electrical  signal.  That  signal  is 
inevitably  noisy.  The  stochastic  models  used  to  characterize  the  noisy  signals  were 
covered  in  Sec.  2.2  of  this  chapter.  Received  LADAR  signals  are  often  used  to  calculate 
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the  range  to  a  target.  This  process,  and  the  limitations  on  range  precision  dne  to 
noise,  were  covered  in  Sec.  2.3.  It  has  been  noted  that  the  constant  of  proportionality 
relating  the  mean  energy  in  the  optical  held  and  the  mean  electrical  signal  (the  gain) 
prodnced  by  an  avalanche  photodiode  can  vary.  This  APD  gain  variation  problem 
was  discnssed  in  Sec.  2.4. 

This  dissertation  covers  several  different  topics  associated  with  the  problem  of 
range  estimation  nsing  a  LADAR  system.  These  items  were  described  in  Sec.  2.5  of 
this  chapter. 
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III.  Poisson-Distributed  LADAR  Signal  Models 

The  subject  of  this  chapter  is  LADAR  signal  modeling  using  the  Poisson  distribution. 
These  models  represent  the  signals  detected  by  LADAR  systems  as  series  of  statisti¬ 
cally  independent  Poisson  random  variables.  The  validity  of  the  Poisson  assumption 
is  discussed  in  detail  in  Sec.  2.2.3.  The  mean  values  of  these  Poisson  random  vari¬ 
ables  are  the  sum  of  a  fixed  bias  level  and  possibly  a  signal  component.  The  bias 
level  represents  the  noise  sources  in  the  system  such  as  background  light  entering 
the  LADAR’s  receiver  aperture  and  dark  current  in  the  present  in  the  photodetector. 
The  signals  are  modeled  as  hnite  duration  (truncated)  parabolas.  The  parabolic  pulse 
shape  model  is  not  used  universally  in  the  literature,  but  is  a  reasonable  approxima¬ 
tion  to  the  true  shape  of  a  laser  pulse.  Various  laser  pulse  shape  models  used  in  the 
literature  are  compared  in  Sec.  2.1.6.  The  parabolic  model  is  used  because  it  is  simple 
and  reasonably  realistic. 

The  LADAR  signal  models  are  parameterized  by  range  to  target,  gain,  bias,  and 
pulse  width.  In  some  cases,  more  than  one  range  or  gain  term  is  present  in  the  signal 
model.  When  multiple  parameters  of  the  same  type  are  used,  it  is  to  simulate  the 
signal  observed  when  the  LADAR  transmits  one  pulse  but  observes  received  signals 
from  more  than  one  target.  It  is  usually  assumed  that  the  received  laser  pulse  width 
is  known.  However,  received  pulse  width  can  be  treated  as  an  unknown  parameter  to 
be  estimated. 

The  LADAR  signal  models  are  used  to  derive  parameter  estimates  and  to  cal¬ 
culate  the  bounds  on  precision  of  those  estimates.  The  use  of  established  parameter 
estimation  techniques  such  as  maximum  likelihood  estimation  is  discussed.  The  signal 
models  are  also  used  to  hnd  the  Cramer-Rao  lower  bounds  for  unbiased  estimates  of 
unknown  parameters. 

Three  different  LADAR  signal  models  are  discussed  in  this  chapter.  The  first 
model  considered  is  the  single  parabolic  pulse  model  with  known  width.  This  model 
was  studied  in  a  previous  publication  [4]  and  is  discussed  in  Sec.  2.2.4.  It  is  revisited 
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in  order  to  present  analytic  expressions  for  the  Cramer-Rao  lower  bounds  on  the 
parameter  estimates  in  that  model. 

The  second  model  in  this  chapter  is  a  generalization  of  the  parabolic  pulse  model 
that  includes  obscuration.  This  obscured  target  model  is  used  to  approximate  the 
observed  LADAR  signal  when  the  target  of  interest  is  partially  occluded  by  clutter. 
In  military  remote  sensing  scenarios,  the  obscurant  is  typically  camouflage  netting  or 
some  sort  of  foliage.  The  obscuration  is  simulated  by  adding  a  second  parabolic  pulse 
to  the  simulated  signal  with  its  own  range  and  gain  values. 

The  hnal  model  is  equivalent  to  the  single  laser  pulse  parabolic  signal  model. 
However,  this  version  treats  the  pulse  width  as  an  unknown  parameter  that  must  be 
estimated  instead  of  a  known  constant.  Estimation  algorithms  for  the  pulse  width 
are  proposed.  The  effect  of  unknown  pulse  width  on  the  Cramer-Rao  lower  bounds  is 
also  discussed. 

3.1  Parabolic  Pulse  LADAR  Signal  Model 

All  three  LADAR  signal  models  in  this  chapter  are  generalizations  of  the  model 
published  by  Cain,  et  al.  [4]  in  'Plash  light  detection  and  ranging  accuracy  limits 
for  returns  from  single  opague  surfaces  via  Cramer-Rao  bounds.”  The  CRLBs  for  the 
parabolic  pulse  model  are  derived  in  that  paper.  However,  they  are  expressed  in  terms 
of  summations  of  various  sampled  functions. 

Because  parabolic  pulse  models  are  used  extensively  in  this  dissertation,  it 
proves  useful  to  derive  analytic  approximations  to  the  summations  from  [4].  These 
analytic  approximations  are  used  to  find  simpler  expressions  for  various  Cramer-Rao 
lower  bounds.  In  this  section,  the  original  parabolic  pulse  model  from  [4]  is  revisited 
in  order  to  derive  analytic  CRLB  expressions. 


35 


In  this  dissertation,  the  round-trip  range  from  Eq.  2.32  is  replaced  by  range  to 


target  R.  Round-trip  range  R^t  is  double  the  one-way  range  to  the  target. 


Rrt  =  2R. 


(3.1) 


The  parabolic  pulse  model  in  terms  of  range  to  target  is 
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The  pulse  is  sampled  in  its  entirety  as  long  as  pw  <  2R/c  <  td  —  Pw 

The  parabolic  pulse  model  from  Eq.  3.2  is  equal  to  Eq.  2.32  with  the  substitu¬ 
tion  from  Eq.  3.1.  An  example  of  the  parabolic  pulse  model  from  Eq.  3.2  is  shown 
in  Fig.  3.1.  The  signal  shapes  shown  in  Figs.  2.3  and  3.1  appear  equal  because  the 
range  to  target  is  the  same  in  both.  In  Fig.  2.3,  the  round-trip  range  is  Rrt  =  10  m 
and  in  Fig.  3.1,  the  range  is  i?  =  5  m. 


3.1.1  Fisher  Information  Matrix.  The  Fisher  information  matrix  is  used  to 
calculate  the  Cramer-Rao  lower  bounds  on  unbiased  estimates  of  signal  parameters. 
The  unknown  parameters  in  Eq.  3.2  are  range,  gain,  and  bias  (denoted  R,  G,  and 
B).  The  FIM  associated  with  the  parabolic  pulse  model  from  Eq.  3.2  is 
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The  FIM  elements  are  calculated  from  the  expectations  of  the  second  derivatives  of  the 
log-likelihood  function  /(R,  G,  B).  The  signal  is  modeled  as  a  sequence  of  independent 
Poisson  random  variables,  each  with  mean  I [tk).  The  log-likelihood  function  from  (see 
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Parabolic  Pulse  LADAR  Signal  Model 


Figure  3.1;  Illustration  of  the  parabolic  pulse  LADAR  signal 
model  from  Eq.  3.2.  This  signal  model  differs  from  the  one 
shown  in  Fig  2.3  because  it  is  defined  in  terms  of  range  to  target 
{R)  rather  than  round-trip  range  {Rrt)- 
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Eqs.  2.35  and  2.36)  for  the  parabolic  pulse  from  Eq.  3.2  is 
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Derivatives  of  the  log-likelihood  function  are  needed  to  calculate  the  FIM.  The 
first  derivatives  of  the  log-likelihood  function  for  the  parabolic  pulse  signal  model  (Eq. 
3.2)  with  respect  to  the  unknown  parameters  are 
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Eq.  3.7  was  evaluated  using  the  fact  that 
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For  the  CRLB  to  exist,  the  expectations  of  the  first  derivatives  of  the  log- 
likelihood  function  must  be  zero.  This  requirement  is  called  the  regularity  condition 
[23].  The  fact  that  the  mean  of  Poisson  random  variable  Ditk)  is  I{tk)  is  used  to  find 
these  expectations. 
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Eq.  3.9  can  be  used  to  show  that  the  first  derivatives  of  the  log-likelihood  (Eqs.  3.5 
-  3.7)  all  have  an  expected  value  of  zero.  Thus,  the  regularity  condition  is  satisfied. 
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The  second  derivatives  needed  for  the  diagonal  elements  of  the  Fisher  informa¬ 


tion  matrix  are 
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The  second  derivatives  needed  for  the  off-diagonal  elements  of  the  FIM  are 
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where  Eq.  3.8  was  used  to  evaluate  Eqs.  3.14  and  3.15. 


The  expectations  of  the  log-likelihood  function  derivatives  must  be  calculated 
to  hnd  the  FIM.  Eq.  3.9  is  used  to  find  the  FIM  element  expectations.  The  FIM 
diagonal  elements  are 
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and 
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The  off-diagonal  FIM  elements  are 
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The  derivatives  of  the  parabolic  pnlse  are  needed  to  evalnate  Eqs.  3.16  -  3.21. 
The  derivative  of  Eq.  3.2  with  respect  to  range  is 
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The  derivative  with  respect  to  gain  is 
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In  these  two  derivative  eqnations  (Eqs.  3.22  and  3.23),  the  effect  of  the  discontinnity 
in  I{tk)  dne  to  the  rectangle  function  has  been  ignored. 

Using  Eq.  3.22,  the  FIM  element  from  Eq.  3.16  can  be  written 
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The  Fisher  information  matrix  element  in  Eq.  3.24  includes  a  summation.  Summa¬ 
tions  such  as  this  one  can  be  approximated  using  integrals.  The  integral  approxima¬ 
tion  to  Eq.  3.24  is 
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This  approximation  assumes  that  the  data  was  sampled  at  rate  fs-  The  integral  in 
Eq.  3.25  is  evaluated  using  the  following  change  of  variable. 
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Using  this  change  of  variable  in  the  integral  from  Eq.  3.25  yields 
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Eq.  2.147.5  of  Table  of  Integrals,  Series,  and  Products  by  Gradshteyn  and  Ryzhik  [14] 
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Letting  m  =  2  in  the  previous  equation  produces 
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1 


(3.30) 


Eq.  2.143.2  of  [14]  is 
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if  |x|  <  1  where  atanh(a;)  is  the  inverse  hyperbolic  tangent  fnnction.  Combining  Eqs. 
3.25,  3.29,  and  3.30  prodnces 
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The  following  variable  a  is  dehned  to  simplify  the  Fisher  information  matrix  elements 
and  later  CRLB  equations. 
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Using  this  variable  allows  Eq.  3.31  to  be  written 
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This  integral  approximation  for  the  FIM  element  was  validated  by  comparing  it  to 
an  actual  summation.  The  results  are  shown  in  Fig.  3.2a.  There  is  good  agreement 
between  the  integral  approximation  and  the  numeric  result. 
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Figure  3.2:  (a)  The  range  second  derivative,  (b)  The  range-gain  derivative,  (c)  The 

range-bias  derivative,  (d)  The  gain  second  derivative,  (e)  The  gain-bias  derivative, 
(f)  The  bias  second  derivative. 
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Using  Eqs.  3.22  and  3.23,  the  FIM  element  from  Eq.  3.19  can  be  approximated 


by  an  integral. 
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The  following  change  of  variable  was  used  in  the  above  equation. 
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Pw 
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Inspection  of  the  integrand  of  Eq.  3.34  shows  that  it  is  an  odd  function  (i.e.,  f{x)  = 
—f{—x)).  Since  the  integrand  is  odd  and  the  range  of  integration  is  centered  around 
the  origin,  the  integral  is  zero.  Therefore,  this  FIM  element  is  approximately  zero. 
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This  approximation  was  compared  to  a  summation  in  Fig.  3.2b.  The  summation 
was  very  close  to  zero,  but  not  exactly  zero.  The  numeric  result  is  not  exactly  zero 
because  the  samples  are  not  aligned  exactly  with  the  pulse  center.  If  they  were,  then 
the  summation  would  be  exactly  zero.  However,  the  magnitude  of  the  summation  is 
very  small  and  can  be  accurately  approximated  by  a  zero. 
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The  FIM  element  from  Eq.  3.20  is  approximately 
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where  the  change  of  variable  from  Eq.  3.35  is  used.  The  integrand  of  Eq.  3.37  is  an 
odd  function  and  the  range  of  integration  is  centered  around  the  origin.  Thus,  the 
integral  is  zero  and  the  FIM  element  is  approximately  zero. 
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The  integrand  being  approximated  is  shown  in  Fig.  3.2c.  In  that  figure,  the  approxi¬ 
mation  was  compared  to  a  summation.  The  summation  result  was  small,  but  not  as 
small  as  in  the  previous  case.  The  discrepancy  is  likely  caused  by  the  sharp  disconti¬ 
nuities  in  the  function  near  the  edges  of  the  parabolic  pulse.  The  other  FIM  element 
with  an  integral  approximation  that  is  approximately  zero  (Eq.  3.36)  is  shown  in  Fig. 
3.2b.  In  that  case,  the  function  did  not  have  sharp  discontinuities  and  the  summation 
was  about  times  smaller  than  the  output  in  this  case.  However,  the  analytic  ap¬ 
proximation  error  is  still  small  compared  to  the  other  terms  in  the  Fisher  information 
matrix.  Thus,  the  approximations  are  accurate  and  the  CRLBs  they  generate  are 
precise. 
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By  using  the  change  of  variable  from  Eq.  3.26,  the  FIM  element  from  Eq.  3.17 


can  be  approximated  by 
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Eq.  3.39  is  the  sum  of  three  integrals.  The  hrst  and  second  integrals  in  this  equation 
can  be  evaluated  using  Eqs.  3.29  and  3.30.  In  order  to  evaluate  the  third  integral  in 
Eq.  3.39,  start  by  letting  m  =  4  in  Eq.  3.28. 
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where  |a;|  <  1  and  Eq.  3.30  is  used.  The  hrst  integral  from  Eq.  3.39  is  evaluated 
using  Eq.  3.30. 
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The  second  integral  from  Eq.  3.39  is  evaluated  using  Eqs.  3.29  and  3.30. 
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The  third  integral  from  Eq.  3.39  is  evaluated  using  Eq.  3.40. 
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Eqs.  3.41,  3.42,  and  3.43  can  be  combined  to  create  an  analytic  approximation  for 
the  Fisher  information  matrix  term  from  Eq.  3.17. 
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Using  the  variable  a,  which  was  dehned  in  Eq.  3.32,  this  can  be  written 
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This  integral  approximation  is  compared  to  a  summation  in  Fig.  3. 2d.  The  numeric 
and  analytic  results  are  exactly  the  same. 

The  FIM  element  from  Eq.  3.21  can  be  approximated  with  an  integral. 
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There  are  two  integrals  in  Eq.  3.46.  These  integrals  can  be  evalnated  using  Eqs.  3.28 
and  3.29.  Therefore,  the  FIM  element  from  Eq.  3.21  is  approximately 
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where  the  variable  a  from  Eq.  3.32  was  used.  It  is  shown  in  Fig.  3.2e  that  this  result 
is  nearly  equal  to  the  summation  over  the  actual  function. 

The  final  FIM  element  to  be  evaluated  (Eq.  3.18)  can  be  approximated  with  an 
integral. 
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This  integral  approximation  differs  from  the  integral  approximations  of  all  the  other 
FIM  terms  because  the  integrand  is  nonzero  over  the  entire  data  collection  region 
(0  <  t  <  td).  This  integral  can  be  computed  by  splitting  it  into  three  different 
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regions. 
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where  the  change  of  variable  from  Eq.  3.26  is  used.  It  can  be  shown  using  Eq.  3.30 
that  Eq.  3.49  is 
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Therefore,  the  FIM  element  from  Eq.  3.18  can  be  approximated  by  the  following 
formula. 


-E 


dH{R,G,B) 


552 


(3.51) 


‘^■Pwfs 


td 


1 

-  + 


—  O  f  i  ‘^Pw 

~  •  2p^B  + 


2p^B  B  ^y/BTG 
a 


atanh 


G 


B  +  G 


B  +  G 


where  Eq.  3.32  was  used.  This  integral  approximation  to  the  FIM  term  is  compared 
to  the  output  of  a  summation  in  3.2f.  The  results  are  nearly  identical. 
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3.1.2  Cramer- Rao  Lower  Bounds.  The  Cramer- Rao  lower  bounds  are  the 
lower  limits  on  the  variance  of  any  unbiased  estimate  of  an  unknown  parameter.  The 
CRLBs  are  obtained  by  inverting  the  Fisher  information  matrix.  The  FIM  for  the 
parabolic  pulse  signal  model  from  Eq.  3.2  is 
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(3.52) 


There  are  matrix  elements  approximately  equal  to  zero  because  Eqs.  3.36  and  3.38 
are  about  zero.  The  Jrr  term  in  the  Fisher  information  matrix  is  Eq.  3.33. 
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where  a  was  dehned  in  Eq.  3.32.  The  Jqg  term  in  the  FIM  is  Eq.  3.45. 
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The  Jgb  term  is  Eq.  3.47. 
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Finally,  the  Jbb  term  is  Eq.  3.51. 
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The  Fisher  information  matrix  (Eq.  3.52)  for  the  parabolic  pulse  LADAR  signal 
model  from  Eq.  3.2  has  a  block  diagonal  form.  The  inverse  of  a  block  diagonal  matrix 
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Since  the  FIM  (Eq.  3.52)  for  the  parabolic  LADAR  pulse  signal  model  from  Eq.  3.2 
is  block  diagonal,  Eq.  3.57  can  be  used  to  hnd  the  FIM  inverse. 


^/•Jrr  0  0 

0  Jss/d  —JcB/d 
0  —JoB/d  Jcc/d 


(3.58) 


where  d  is  the  determinant  of  the  2x2  matrix  in  the  lower  right  corner  of  the  3x3 
FIM  matrix  in  Eq.  3.52.  The  determinant  of  the  2x2  matrix  is 
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Eqs.  3.45,  3.47,  and  3.51  are  used  to  calculate  this  determinant. 
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The  CRLBs  for  estimates  of  range,  gain,  and  bias  are  the  diagonal  elements  of 
the  inverse  of  the  Fisher  information  matrix.  The  CRLB  for  range  estimates  is  the 
reciprocal  of  Jur.  This  variable  is  used  in  Eq.  3.53  to  calculate  the  bound. 
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The  CRLB  for  estimates  of  the  gain  is  Jbb  (Eq.  3.56)  divided  by  the  determinant  d 
from  Eq.  3.60. 
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Table  3.1;  Parabolic  Pulse  Signal 
Model  CRLBs. 


Parameter 

Cramer-Rao  Lower  Bounds 

Range 
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The  bias  estimate  CRLB  is  Jgg  (Eq-  3.54)  divided  by  the  determinant  d  from  Eq. 
3.60. 
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The  following  term  appears  in  the  gain  and  bias  CRLBs  (Eqs.  3.62  and  3.63)  and  can 
be  expressed  in  terms  of  the  hyperbolic  arctangent  function. 
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The  Cramer-Rao  lower  bounds  for  the  parabolic  pulse  LADAR  signal  model  from  Eq. 
3.2  are  summarized  in  Table  3.1. 


3.1.3  Maximum  Likelihood  Estimation  of  Parameters.  Maximum  likelihood 
estimation  is  one  of  the  most  popular  parameter  estimation  techniques.  The  MLEs 
for  the  range,  gain,  and  bias  in  the  parabolic  pulse  signal  model  from  Eq.  3.2  are 
calculated  from  the  joint  PDF  for  the  samples  (Eq.  2.33). 

The  likelihood  function  is  maximized  with  respect  to  the  unknown  values.  Max¬ 
imization  of  the  likelihood  function  is  equivalent  to  maximization  of  the  log-likelihood 
function  from  Eq.  3.4.  Since  the  log-likelihood  function  is  easier  to  work  with  with 
respect  to  range,  it  is  used  to  find  the  MLEs.  The  derivative  of  the  log-likelihood 
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function  with  respect  to  range  is  (see  Eqs.  3.5  and  3.22) 
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The  second  term  has  odd  symmetry  about  the  center  of  the  rectangle  function.  There¬ 
fore,  it  goes  to  zero.  The  MLE  is  found  by  setting  Eq.  3.65  to  zero.  The  range  MLE 
is  the  value  of  R  that  satishes  the  following  equation. 
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Evaluation  of  Eq.  3.66  requires  knowledge  of  the  gain  and  bias.  In  real  LADAR  data, 
the  gain  and  bias  are  usually  unknowns.  Thus,  the  range,  gain,  and  bias  usually  have 
to  be  estimated  simultaneously  in  order  to  calculate  the  range  MLE. 

The  derivative  of  the  log-likelihood  function  with  respect  to  gain  is  (see  Eqs. 
3.6  and  3.23) 
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By  using  the  change  of  variable  from  Eq.  3.35,  the  second  summation  in  Eq.  3.67 
can  be  written 
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Therefore,  the  gain  MLE  is  the  value  of  G  satisfying 
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The  derivative  of  the  log-likelihood  function  with  respect  to  bias  is 
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Therefore,  the  MLE  of  the  bias  is  the  value  B  that  satishes 
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All  three  of  the  MLEs  derived  (Eqs.  3.66,  3.69,  and  3.71)  require  knowledge  of 
the  other  unknown  parameters  in  the  model.  In  practice  all  three  MLEs  must  be  solved 
simultaneously.  Additionally,  there  is  no  analytic  solution  for  any  of  these  MLEs.  All 
must  be  solved  using  iterative  numeric  methods.  Iterative  numeric  estimation  of 
multiple  parameters  is  often  complicated  and  can  be  time  consuming.  Therefore,  in 
Sec.  3.1.4  alternate  methods  of  estimation  are  considered. 
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3.1.4  Range  Estimation  Methods.  It  is  noted  in  Sec.  3.1.3  that  the  maximum 
likelihood  estimates  of  range,  gain,  and  bias  in  the  parabolic  pulse  model  are  difficult 
to  calculate.  Calculation  is  complicated  because  the  MLEs  must  be  evaluated  using 
iterative  numeric  methods  implemented  simultaneously  on  all  three  parameters.  This 
complexity  leads  to  consideration  of  simpler  parameter  estimates. 

It  is  proposed  that  range  estimation  can  be  done  using  a  combination  of  matched 
hltering  and  peak  fitting.  The  problem  of  range  estimation  for  a  signal  that  has 
been  mixed  with  zero-mean  white  Gaussian  noise  was  considered  in  Fundamentals 
of  Statistical  Signal  Processing:  Estimation  Theory  by  Kay  [23] .  It  was  determined 
that  for  a  signal  mixed  with  zero-mean  Gaussian  noise,  that  the  range  MLE  can  be 
calculated  by  filtering  the  noisy  signal  with  a  noiseless  copy  of  the  signal.  The  range 
MLE  is  the  range  at  which  the  hltered  output  is  at  its  maximum. 

Range  estimation  in  LADAR  systems  differs  from  the  zero-mean  Gaussian  noise 
case  because  the  received  signal  is  Poisson-distributed  and  never  has  nonnegative 
values.  Also,  while  it  is  assumed  that  each  sample  in  the  received  signal  is  statistically 
independent  from  every  other  sample,  the  variance  is  changes  from  sample  to  sample. 
In  fact,  the  variance  is  equal  to  the  mean  value  for  Poisson  noise. 

The  range  estimate  is  based  on  the  Gaussian  noise  MLE  from  [23] .  The  estimate 
is  the  range  at  which  the  output  of  the  data  d{t)^  convolved  with  hlter  h{t),  is  at  its 
maximum.  The  range  estimate  R  can  be  written 

R  =  argmax/i(t)  *  d(t).  (3.72) 

R 

Because  the  Poisson-distributed  noise  is  not  white,  the  matched  filter  is  not 
necessarily  the  optimal  hlter  for  range  estimation.  Thus,  alternative  Liters  are  tested 
in  the  range  estimation  process.  Four  different  filters  were  tested  in  Eq.  3.72.  In  the 
hrst  case,  the  filter  is  equal  to  the  shape  of  the  parabolic  pulse  Pp(t).  That  is, 

hiit)=Ppit)  (3.73) 
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where  the  pulse  is 
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This  is  a  true  matched  filter  where  hi{t)  is  equal  to  the  noiseless  pulse.  The  next  filter 
tested  was  the  square  root  of  the  noiseless  pulse.  The  reason  for  trying  this  variation 
on  the  true  matched  filter  is  that  in  the  Poisson  noise  case,  the  samples  with  the 
highest  mean  also  have  the  greatest  variance.  Using  the  square  root  of  the  noiseless 
pulse  will  put  slightly  less  weight  on  the  highest  variance  samples  than  the  standard 
matched  filter  does.  The  filter  is, 


Hi)  =  |Pp(i)l''".  (3.75) 

The  third  filter  tested  uses  the  true  mean  value  of  the  received  signal.  Because  the 
noise  (the  bias)  is  not  zero-mean,  the  bias  value  is  added  to  the  filter.  This  filter  is 

h^it)  =  Gpp{t)  +  B.  (3.76) 

The  final  filter  tested  was  the  square  root  of  the  true  mean  value.  The  motivation 
for  using  this  filter  is  the  same  as  for  filter  h2{t)  (Eq.  3.75).  That  is,  the  square  root 
operation  puts  slightly  less  weight  on  the  highest  variance  values  in  the  received  signal 
than  the  filter  without  the  square  root  (Eq.  3.76)  does.  This  filter  is. 

Kit)  =  [Gppit)  +  BY/\  (3.77) 

There  is  a  fundamental  difference  between  the  first  two  and  the  last  two  filters. 
The  first  two  filters  (Eqs.  3.73  and  3.75)  do  no  require  any  knowledge  of  the  unknown 
parameters  in  the  system.  The  second  two  filters  (Eqs.  3.76  and  3.77)  assumed  exact 
knowledge  of  the  gain  and  bias  in  the  system.  In  practice,  the  gain  and  bias  may  not  be 
known  exactly  and  may  need  to  be  replaced  with  estimates  in  order  to  implement  the 
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second  two  filters.  Filters  that  assume  knowledge  of  system  parameters  that  could  be 
unknowns  are  sometimes  called  clairvoyant  hlters  [10].  Because  the  second  two  filters 
assume  knowledge  of  the  gain  and  bias  in  the  system,  they  are  called  clairvoyant. 

The  four  filters  were  tested  using  Monte  Carlo  simulations.  The  results  are 
shown  in  Fig.  3.3.  Simulated  Poisson  noise  signals  were  generated.  A  simulated 
ten  nanosecond  pulse  =  10  ns)  was  centered  within  a  one  hundred  nanosecond 
{td  =  100  ns)  sampling  interval.  This  signal  was  sampled  at  a  rate  of  one  Gigahertz 
(/s  =  1  GHz).  Therefore,  there  were  one  hundred  samples  in  each  simulated  received 
signal  {K  =  100).  The  bias  in  the  signal  was  fixed  at  five  {B  =  5)  and  the  gain 
was  varied  from  one-tenth  to  one  thousand  (0.1  <  G  <  1000).  The  Monte  Garlo 
simulation  was  run  many  times  in  order  to  obtain  multiple  range  estimates  for  each 
gain  value  simulated.  The  precision  of  the  range  estimate  was  measured  by  calculating 
the  sample  mean  over  the  set  of  range  estimates  for  each  of  the  four  hlters. 

The  simulated  received  signals  were  processed  using  each  of  the  four  hlters  (Eqs. 
3.73  -  3.77).  The  simulated  signal  was  convolved  with  each  of  the  hlters,  then  the 
edges  of  the  output  were  trimmed  so  that  the  width  of  the  sampling  interval  was 
not  increased.  The  range  estimate  was  obtained  by  hnding  the  peak  sample  and 
then  doing  a  three-point  parabolic  ht  to  the  peak  output  point  and  the  two  adjacent 
points.  The  exact  parabolic  ht  algorithm  used  was  the  one  presented  in  Eqs.  1  and  2 
of  [20].  The  parabolic  ht  was  used  to  get  precision  better  than  a  bin  width  from  the 
range  estimate.  If  the  range  estimate  had  been  set  to  the  center  of  the  peak  bin,  then 
the  precision  could  never  be  better  than  the  standard  deviation  of  a  uniform  random 
variable  over  one  bin  width. 

The  sampling  duration  is  td-  The  measured  round-trip  time  trt  is  related  to  the 
round-trip  range  Rrt  by 

Rrt  =  UtC  (3.78) 

where  c  is  the  speed  of  light.  The  variance  of  a  uniform  random  variable  is  the  width 
of  that  random  variable’s  support  squared  divided  by  twelve.  Therefore,  the  variance 
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of  a  range  estimate  that  is  uniformly  distributed  over  a  sampling  duration  of  td  is 


48  ■ 


(3.79) 


Inspection  of  the  results  shown  in  Fig.  3.3  shows  that  none  of  the  filters  tested 
achieved  the  Cramer-Rao  lower  bound.  However,  for  gain  values  between  about  three 
and  one  thousand,  the  best  precision  was  at  worst  about  double  the  CRLB.  As  the 
gain  values  dropped  below  the  level  of  the  bias,  the  range  precision  decreased  rapidly. 
For  gains  around  one  and  two,  the  range  estimate  standard  deviation  is  nearly  ten 
times  the  CRLB.  When  the  gain  is  low  enough  (less  than  about  one),  the  range 
estimate  precision  is  equal  to  the  precision  of  a  uniform  random  estimate  over  the 
search  interval  (see  Eq.  3.79).  A  uniform  random  error  is  effectively  the  worst  case 
scenario.  Therefore,  the  hltering  operation  is  providing  no  benefit  for  gain  values 
below  one. 

The  Ziv-Zakai  lower  bound  (ZZLB)  is  sometimes  used  to  characterize  limits 
on  parameter  precision.  The  ZZLB  was  originally  derived  in  [57]  and  was  improved 
in  [6].  In  [24],  it  was  applied  to  estimation  of  pulse  time-of-flight.  The  ZZLB  has  the 
advantage  of  being  valid  in  the  low  signal-to-noise  (SNR)  regime.  However,  it  is  more 
complicated  to  evaluate  than  the  CRLB  [24].  Also,  the  ZZLB  is  not  as  tight  as  the 
CRLB  at  high  SNRs  [6].  Thus,  it  was  not  used  in  this  dissertation. 

The  estimate  variances  in  Fig.  3.3  drop  below  the  variance  predicted  by  the 
CRLB  because  at  very  low  gain  values,  the  are  biased.  Unbiased  parameter  esti¬ 
mates  can  not  have  variances  lower  than  what  is  predicted  by  their  Cramer-Rao  lower 
bounds.  However,  it  is  possible  for  a  biased  estimate  to  have  a  variance  lower  than 
the  CRLB. 

The  best  filter  overall  is  the  square  root  hlter  (Eq.  3.75).  However,  the  standard 
matched  hlter  (Eq.  3.73)  did  perform  slightly  better  than  the  square  root  hlter  for  a 
small  range  of  gain  values.  The  standard  matched  hlter  worked  the  best  for  gains  from 
about  three  to  ten.  The  clairvoyant  hlters  (the  hlters  that  assumed  exact  knowledge  of 


59 


Range  Est.  (p^=  10ns  B=5  f^=1GHz  t^=0.1)a.s) 


Figure  3.3:  Illustration  shows  range  precision  produced  by 

four  different  filters.  The  Cramer-Rao  lower  bound  is  also 
shown. 

the  gain  and  bias)  did  not  work  better  than  the  standard  matched  hlter  and  the  square 
root  matched  hlter.  In  general,  the  precision  of  parameter  estimation  improves  when 
more  of  the  other  parameters  are  known.  But  in  this  case,  knowledge  of  the  gain  and 
bias  did  not  help.  However,  it  is  possible  that  a  different  hlter  that  incorporated  the 
gain  and  bias  values  in  a  diherent  way  could  achieve  performance  closer  to  the  Cramer- 
Rao  lower  bound.  The  clairvoyant  hlters  performed  worse  than  their  counterparts. 
This  performance  convergence  was  expected  since  as  G/B  approaches  inhnity,  the 
hlters  themselves  converge  to  the  non-clairvoyant  versions. 

Several  more  Monte  Carlo  simulations  of  range  estimation  are  shown  in  Fig.  3.4. 
Since  the  results  in  Fig.  3.3  showed  that  the  square  root  hlter  had  the  best  overall 
performance,  it  was  the  only  hlter  used  in  these  simulations.  In  these  simulations, 
the  sampling  duration  was  td  =  100  ns  and  the  sampling  rate  was  /s  =  1  GHz. 
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Range  Est.  (B=5fg=1GHzt^=0.1|is) 


Range  Est.  (p^=10ns fg=1GHzt^=0.1)is) 


Figure  3.4:  (a)  Range  estimation  with  bias  B  =  5. 

(b)  Range  estimation  with  pulse  width  =  10  ns. 

In  Fig.  3.4a,  the  bias  is  R  =  5  and  the  plot  is  parametric  in  pulse  width.  The 
narrowest  pulse  {p^  =  5  ns)  had  the  best  performance  over  most  of  the  gain  values. 
With  all  else  being  equal,  short  pulses  are  better  than  long  pulses  for  range  estimation 
since  the  received  signal  energy  is  focused  within  a  shorter  duration.  When  a  long 
pulse  is  used,  the  energy  is  spread  over  a  longer  duration  and  it  is  more  difficult  to 
obtain  an  precise  measurement  of  the  peak  value  (the  range  estimate).  However,  the 
longest  pulse  =  20  ns)  shows  performance  closest  to  the  CRLB  and  showed  the 
best  performance  at  certain  low  gain  values  (below  about  five).  It  is  possible  that  the 
longer  pulse  is  closer  to  the  CRLB  than  the  shorter  pulses  because  the  shorter  pulses 
have  more  energy  at  higher  frequencies  and  are  not  being  sampled  quickly  enough. 
The  Nyquist  frequency  (/„  =  /s/2)  was  500  MHz  in  these  simulations.  The  5  ns  pulse 
has  four  times  the  bandwidth  of  the  20  ns  pulse,  so  it  is  possible  that  the  20  ns  pulse 
is  being  sampled  at  an  adequate  rate  while  the  5  ns  pulse  is  experiencing  aliasing 
when  sampled. 

The  bandwidth  of  a  parabolic  pulse  was  derived  in  Appendix  B.  It  is  shown  in 
Fig.  B.2  that  only  about  two  parts  in  ten  thousand  of  the  total  energy  is  outside  the 
Nyquist  frequency  for  the  5  ns  pulse.  The  energy  outside  the  Nyquist  frequency  is 
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aliased  and  will  interfere  with  parameter  estimates.  The  fraction  of  parabolic  pulse 
energy  aliased  for  the  10  and  20  ns  pulses  is  nearly  one-tenth,  and  one  one-hundredth, 
lower  than  for  the  5  ns  pulse,  respectively.  The  lower  amount  of  aliased  energy  is  a 
possible  explanation  for  why  the  longer  pulses  were  closer  to  the  range  estimate  CRLB 
in  Fig.  3.4a. 

In  Fig.  3.4b,  the  pulse  width  is  =  10  ns  and  the  plot  is  parametric  in  bias. 
It  is  not  surprising  that  the  best  range  estimate  precision  was  observed  when  the  bias 
was  lowest  {B  =  1).  For  bias  B  =  1,  the  range  estimate  was  nearly  at  the  CRLB  for 
gain  values  over  about  three.  For  the  highest  bias  case  {B  =  10),  the  precision  was 
signihcantly  over  the  CRLB  until  the  gain  was  greater  than  about  six.  Once  the  gain 
was  over  about  one  hundred,  the  range  precision  was  roughly  equal  for  all  three  bias 
values  simulated  because  at  those  levels,  the  bias  was  much  smaller  than  the  gain  in 
each  case. 


3.1.5  Gain  and  Bias  Estimation.  In  Sec.  3.1.3,  the  maximum  likelihood 
estimates  for  gain  (Eq.  3.69)  and  bias  (Eq.  3.71)  are  derived.  However,  the  MLEs  for 
both  of  these  parameters  are  difficult  to  evaluate.  In  this  section,  simple  clairvoyant 
estimates  of  gain  and  bias  are  proposed  and  tested.  The  estimates  are  clairvoyant 
because  they  assume  knowledge  of  the  range  to  the  target.  These  clairvoyant  estimates 
are  not  useful  with  real  data  since  they  assume  knowledge  of  a  parameter  that  is 
unknown  (the  range).  The  purpose  of  these  estimates  is  to  validate  the  Cramer- Rao 
lower  bound  formulas  derived  in  Sec.  3.1.2. 

The  bias  estimate  used  throws  out  the  points  in  the  received  data  that  contain 
any  part  of  the  parabolic  pulse.  The  estimator  then  averages  the  remaining  points  to 
obtain  the  bias  estimate.  Mathematically,  the  bias  estimate  is 
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Bias  Estimation  (B=5  f^=1GHz  t^=0.1|is) 


Bias  Estimation  (p^=20ns  f^=1GHz  t^=0.1ps) 


Figure  3.5:  (a)  Bias  estimation  with  bias  B  =  5. 

(b)  Bias  estimation  with  pulse  width  =  10  ns. 

The  same  estimate  was  used  in  [4].  Eq.  3.80  is  equal  to  Eq.  3  from  that  reference. 
This  estimator  is  simply  an  average  of  the  data  points  that  are  separate  from  the 
pulse.  The  rect(a;)  term  in  the  numerator  and  denominator  zero  the  points  that  are 
part  of  the  parabolic  pulse  in  the  received  signal. 

Monte  Carlo  simulations  of  gain  estimation  are  shown  in  Fig.  3.6.  As  in  the 
range  estimation  simulations  (Figs.  3.3,  3.4),  the  sampling  duration  is  td  =  100  ns 
and  the  sampling  rate  is  /^  =  1  GHz.  In  Fig.  3.6a,  the  bias  is  fixed  at  i?  =  5  and  the 
plot  is  parametric  in  pulse  width.  In  Fig.  3.6b,  the  pulse  width  is  fixed  at  p^  =  10 
ns  and  the  plot  is  parametric  in  bias.  In  both  Fig.  3.6a  and  3.6b,  the  precision  of 
the  bias  estimate  (Eq.  3.80)  meets  the  CRLB  unless  the  pulse  duration  is  long  and 
the  gain  is  low.  If  the  pulse  duration  is  short,  then  very  little  data  is  omitted  by  the 
rect(a;)  function  in  the  estimator.  If  the  gain  is  high,  then  the  data  in  the  pulse  is  of 
little  use  for  bias  estimation  regardless  of  the  pulse  duration. 


63 


A  relatively  simple  gain  estimate  can  be  derived  by  first  examining  the  expec¬ 


tation  of  the  sum  of  the  data  points  from  the  pulse. 
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where  the  change  of  variable  from  Eq.  3.35  was  used.  Therefore,  the  following  equa¬ 
tion  is  an  unbiased  gain  estimate. 
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where  the  bias  estimate  from  Eq.  3.80  was  used.  The  rect(a:)  function  in  this  estimate 
omits  the  points  that  are  outside  the  parabolic  pulse.  This  estimator  was  tested  using 
Monte  Carlo  simulations  and  the  results  are  shown  in  Fig.  3.6.  In  Fig.  3.6a,  the 
bias  is  hxed  at  5  =  5  and  the  plot  is  parametric  in  pulse  width.  The  precision  of 
the  estimate  improves  as  the  pulse  width  increases  since  an  increase  in  pulse  width 
provides  more  samples  to  use  in  the  estimator.  In  Fig.  3.6b,  the  pulse  width  is  hxed 
at  Pw  =  10  ns  and  the  bias  is  varied.  The  precision  of  the  gain  estimator  improves  as 
the  bias  drops  since  a  decrease  in  bias  decreases  the  variance  in  the  samples  used  in 
the  estimator. 
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Gain  Estimation  (B=5  fg=1GHz  t^=0.1|j.s) 


Gain  Estimation  (p^=10ns fg=1GHzt^=0.1iis) 


(a) 


(b) 


Figure  3.6:  (a)  Gain  estimation  with  bias  B  =  5. 

(b)  Gain  estimation  with  pulse  width  =  10  ns. 


3.2  Obscured  Target  LADAR  Signal  Model 

The  LADAR  signal  model  published  in  [4]  is  a  representation  of  a  received  signal 
from  a  target  in  the  open.  However,  imaging  LADARs  are  often  used  to  create  images 
of  targets  that  are  behind  vegetation  or  beneath  tree  canopy  [26,27].  In  this  section, 
a  new  signal  model  is  developed  that  includes  the  target  and  the  obscuration.  This 
new  model  is  used  to  derive  estimators  for  parameters  including  the  target  range.  The 
model  is  also  used  to  derive  GRLBs  for  the  signal  parameters. 


3.2.1  LADAR  Signal  Model  with  Obscuration.  In  order  to  develop  efficient 
obscured  target  detection  algorithms,  a  signal  model  is  needed  that  includes  the  ob¬ 
scuration  and  the  target.  The  model  in  [4]  can  be  modihed  to  include  two  returns, 
one  from  the  obscuration  and  one  from  the  target  behind  the  obscurant.  There  are 
different  range  {R)  and  gain  (G)  values  for  both  of  the  returns.  As  with  the  single 
return  model,  the  bias  is  assumed  to  be  at  a  constant  level  over  all  samples.  The 
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following  modification  of  the  parabolic  pulse  model  from  Eq.  3.2  is  proposed. 
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rect 
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‘^Pw 


+  B. 


(3.83) 


This  equation  is  similar  to  the  single  return  parabolic  pulse  model  from  Eq.  3.2. 
However,  it  has  separate  gain  and  range  values  for  the  obscuration  (Gi  and  Ri)  and 
for  the  target  {G2  and  R2).  It  is  assumed  that  Ri  is  the  range  to  the  object  closer  to 
the  sensor  and  that  R2  is  the  range  to  the  farther  object  (i?i  <  i?2)- 

It  is  assumed  that  the  two  pulses  do  not  overlap  and  that  both  pulses  are 
completely  sampled.  In  order  for  the  pulses  to  not  overlap,  it  is  required  that  R2—R1  > 
PyjC.  The  requirements  that  must  be  satished  for  the  pulses  to  be  completely  sampled 
are  2Ri/c  >  p^  and  2R2/C  <  t^  —  Pw-  At  the  speed  of  light,  every  nanosecond  of 
pulse  width  corresponds  to  about  one  foot  of  distance.  Pulses  that  are  only  a  few 
nanoseconds  are  able  to  resolve  multiple  targets  that  are  separated  by  a  few  feet. 

The  obscured  target  signal  model  from  Eq.  3.83  is  illustrated  in  Fig.  3.7.  As  in 
the  original  parabolic  pulse  signal  model  from  Eq.  3.2,  is  the  three-quarter  width 
of  the  maximum  pulse  height.  The  entire  pulse  width  is  2pyj.  The  bias  illustrated 
in  Fig.  3.7  is  one-half  {B  =  0.5).  There  are  two  signals.  One  represents  the  signal 
received  from  the  obscuration.  The  other  signal  is  from  the  target.  The  obscuration 
range  and  gain  illustrated  in  Fig.  3.7  are  =  5  m  and  Gi  =  2.  The  target  range 
and  gain  are  =  10  m  and  G2  =  0.25. 
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Obscured  Target  LADAR  Signal  Model 


Figure  3.7:  Illustration  of  the  obscured  target  LADAR  signal 
model  with  parabolic  pulses.  This  model  can  be  used  for  CRLB 
derivations  and  to  calculate  parameter  estimates. 
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3.2.2  Fisher  Information  Matrix  Elements.  The  Fisher  information  matrix 
for  the  variable  pulse  width  model  from  Eq.  3.83  is 
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(3.84) 


The  log-likelihood  derivatives  with  respect  to  Ri  and  Gi  are  equal  to  zero  outside 
the  range  \t  —  2Ri/c\  <  (see  Eqs.  3.22  and  3.23).  Therefore,  the  second  derivatives 
from  the  FIM  (Eq.  3.84)  are  zero  for  any  pair  of  parameters  that  characterizes  the 
obscuration  and  the  target.  Thus,  there  are  several  FIM  terms  with  expected  value 
zero. 
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The  range  and  gain  second  derivatives  for  the  signal  from  the  obscuration  and  target 
will  also  be  about  zero  as  they  were  for  the  single  parabolic  pulse  case  (see  Eq.  3.36). 
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The  range  and  bias  second  derivatives  are  also  about  zero  (see  Eq.  3.38). 
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Therefore,  the  only  nonzero  off-diagonal  terms  in  the  Fisher  information  matrix  are 
the  second  derivative  pairs  corresponding  to  the  two  gain  variables  and  the  bias. 


The  diagonal  elements  for  range  and  gain  are  identical  to  the  ones  derived  for 
the  parabolic  pulse  model  in  Sec.  3.1.  The  obscured  target  model  FIM  (Eq.  3.84)  can 
be  populated  using  Eqs.  3.33  and  3.45.  The  range  terms  are 
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The  gain  terms  are 
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The  approximation  to  the  FIM  bias  element  is  not  the  same  as  the  one  in  the 
parabolic  pulse  model  because  of  the  multiple  pulses  present  in  the  received  signal. 
However,  it  is  evaluated  in  the  same  way  as  was  done  in  Eqs.  3.48  -  3.51.  Combining 
Eqs.  3.18  and  3.83  yields 
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Through  the  change  of  variable  from  Eq.  3.26,  this  equation  can  be  written 
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The  integrals  in  the  previous  equation  are  equivalent  to  the  integral  that  appears  in 
Eq.  3.50.  Using  the  result  from  that  equation,  the  FIM  term  can  be  written 
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This  analytic  approximation  is  compared  to  the  output  of  a  numeric  integration  in 
Fig.  3.8b  and  the  results  are  nearly  identical. 


The  only  nonzero  off-diagonal  FIM  term  is  the  one  that  corresponds  to  the  gain 
and  the  bias.  This  term  is  identical  to  the  one  from  the  parabolic  pulse  model.  It  is 
(see  Eq.  3.47) 
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3.2.3  Cramer- Rao  Lower  Bounds.  The  Cramer- Rao  lower  bounds  are  the 
lower  limits  on  the  variance  of  any  unbiased  estimate  of  an  unknown  parameter.  The 
CRLBs  are  obtained  by  inverting  the  Fisher  information  matrix.  The  FIM  for  the 
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The  d  function 


parabolic  pulse  signal  model  from  Eq.  3.83  is 
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Several  matrix  elements  approximately  equal  to  zero  because  the  FIM  elements  shown 
in  Eqs.  3.85,  3.86,  and  3.87  are  about  zero.  The  nonzero  elements  of  the  approximate 
FIM  from  Eq.  3.95  are 
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and 
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where  the  variable  was  defined  in  Eq.  3.89. 


(3.99) 


The  Cramer-Rao  lower  bounds  for  the  unknown  parameters  in  the  obscured 
target  signal  model  (Eq.  3.83)  are  the  diagonal  elements  of  the  inverse  of  the  matrix 
in  Eq.  3.95.  The  FIM  is  block  diagonal,  so  it  is  inverted  using  Eq.  3.57.  The 
range  CRLBs  are  the  reciprocals  of  the  corresponding  FIM  elements.  The  CRLB  for 
estimates  of  range  is 
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The  CRLBs  for  gain  and  bias  are  more  complicated  than  the  range  CRLB  since 
there  are  off-diagonal  elements  that  complicate  the  matrix  inverse.  The  gain  CRLBs 
are 
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Finally,  the  bias  CRLB  is 
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The  CRLBs  for  the  obscured  target  LADAR  signal  model  from  Eq.  3.83  are 
summarized  in  Tables  3.2  and  3.3.  The  actual  CRLBs  are  shown  in  Table  3.3.  The 
Fisher  information  matrix  elements  needed  to  calculate  those  CRLBs  are  shown  in 
Table  3.2. 


3.2.4  Parameter  Estimation.  Estimation  of  the  five  parameters  in  the  ob¬ 
scured  target  model  from  Eq.  3.83  is  difficult.  It  is  straightforward  to  express  the 
likelihood  function  associated  with  Eq.  3.83.  However,  there  are  no  simple  analytic 
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Table  3.2:  Obscured  Target  Model 

Fisher  Information  Matrix  Elements. 


Table  3.3:  Obscured  Target  Model 

Cramer-Rao  Lower  Bounds. 
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solutions  that  simultaneously  maximize  the  unknowns  in  that  likelihood  function.  A 
numeric  search  algorithm  such  as  a  gradient  descent  based  search  would  be  necessary 
to  calculate  the  MLEs. 

It  is  desirable  to  test  the  CRLBs  from  Sec.  3.2.3  using  Monte  Carlo  methods. 
Instead  of  relying  on  numeric  methods,  simple  suboptimal  estimators  were  used  on 
simulated  data  following  the  obscured  target  model.  These  estimators  are  similar  to 
the  ones  presented  in  Sec.  3.1.4  and  Sec.  3.1.5,  but  are  more  complicated  because  the 
received  signal  now  contains  two  different  signals. 

Four  different  range  estimator  hlters  were  tested  using  Monte  Carlo  simulations 
in  Sec.  3.1.4.  The  filter  with  the  best  overall  performance  was  h2{t)  (see  Eq.  3.75). 
This  hlter  is  used  to  estimate  the  range  to  the  largest  gain  signal.  The  estimate  is 


Riarge  =  argmaxh2(t)  *  d(t) 
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(3.104) 


where  hlter  h2{t)  is  the  square  root  of  the  parabolic  pulse  from  the  obscured  target 
model  (Eq.  3.83).  This  hlter  is 
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The  estimate  above  is  the  range  to  the  peak  gain  signal.  However,  there  are 
two  signals  in  the  obscured  target  model.  The  range  to  the  other  signal  is  calculated 
using  the  same  hlter  as  in  Eq.  3.104.  However,  the  large  gain  pulse  centered  at  range 
Riarge  IS  removed  from  the  data  before  the  hlter  is  employed.  The  range  estimate  for 
the  smaller  gain  pulse  can  be  expressed 
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The  obscured  target  signal  illustrated  in  Fig.  3.7  shows  gain  being  larger  for 
the  closest  pulse  (Gi  >  ^2)-  However,  it  was  not  assumed  that  that  the  closest 
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pulse’s  gain  is  always  larger.  Therefore,  the  range  estimate  to  the  closest  pulse  is  the 
minimum  of  the  estimates  from  Eq.  3.104  and  Eq.  3.106. 


Ri  =  min 


(3.107) 


The  range  estimate  to  the  pulse  that  is  farther  away  is 


(3.108) 


As  in  Sec.  3.1.5,  simple  clairvoyant  gain  and  bias  estimates  are  used  in  order 
to  validate  the  CRLBs  derived  in  Sec.  3.2.3.  The  estimates  are  clairvoyant  because 
they  assume  exact  knowledge  of  the  target  ranges  Ri  and  i?2-  The  bias  estimate  is  a 
generalization  of  the  one  shown  in  Eq.  3.80.  It  is 


(3.109) 


This  bias  estimate  is  similar  to  the  one  in  Eq.  3.80,  but  in  this  case,  there  are  two 
rectangle  functions.  Two  rectangle  functions  are  needed  to  remove  the  two  signals  in 
the  obscured  target  model. 

The  clairvoyant  gain  estimates  are  based  on  the  one  shown  in  Eq.  3.82.  Like 
the  bias  estimate  from  Eq.  3.109,  they  assume  exact  knowledge  of  the  target  ranges 
Ri  and  i?2-  They  are 
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and 
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3.2.5  Simulations.  Range  estimation  for  an  obscnred  target  was  simnlated 
using  Monte  Carlo  methods.  Many  realizations  of  Poisson-distributed  signals  with 
mean  values  given  by  Eq.  3.83  were  generated.  Using  the  randomly  generated  data, 
the  ranges  to  the  near  and  far  targets  were  estimated  using  Eqs.  3.107  and  3.108. 
These  range  estimate  equations  do  not  require  knowledge  of  the  gain  or  bias  in  the 
data.  Thus,  they  are  not  clairvoyant  estimators.  The  results  of  the  simulation  are 
shown  in  Fig.  3.9. 

This  obscured  target  range  estimation  simulation  is  similar  to  the  range  esti¬ 
mation  simulation  from  the  previous  section  (Figs.  3.3  and  3.4).  However,  there 
are  some  differences  between  the  range  estimation  simulations  in  the  different  sec¬ 
tions.  As  the  gain  drops,  the  range  standard  deviation  in  the  simulation  shown  in 
Fig.  3.9  settles  to  value  that  is  less  than  the  uniformly-distributed  error  from  Eq. 
3.79.  This  uniformly-distributed  error  is  not  achieved  because  the  range  estimation 
equations  (Eqs.  3.107  and  3.108)  bias  the  estimates  toward  the  beginning  and  end  of 
the  interval,  respectively.  At  high  gain  values,  the  CRLB  is  nearly  achieved. 

Monte  Carlo  simulations  of  gain  estimation  are  shown  in  Fig.  3.10.  In  Fig. 
3.10a,  the  gain  of  the  hrst  pulse  (Gi)  is  estimated  while  the  gain  of  the  second  (G2) 
is  varied.  These  simulations,  and  the  CRLBs  illustrated  in  the  hgure,  show  that  the 
precision  of  estimates  of  one  gain  value  are  independent  of  the  other  gain  value.  In 
Fig.  3.10b,  the  gain  of  the  second  pulse,  G2,  is  estimated  for  different  values  of  that 
gain.  As  G2  increases,  the  variance  of  its  estimate  also  increases.  This  occurs  since 
the  data  is  Poisson-distributed  and  the  variance  increases  with  increases  in  gain  or 
bias.  However,  given  a  particular  value  of  R,  the  CRLB  only  drops  so  far  because  the 
variance  of  the  bias  remains  to  interfere  with  the  precision  of  the  estimate. 

Monte  Carlo  simulations  of  bias  estimation  are  illustrated  in  Fig.  3.11.  In  Fig. 
3.11a,  the  gain  of  the  second  pulse,  G2,  is  varied  and  the  simulations  are  run  for  three 
different  bias  values  {B  =1,  3,  and  10).  The  bias  estimate  variance  is  dominated  by 
the  actual  bias  level.  The  CRLB  does  show  a  slight  dependence  on  the  gain  level. 
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standard  Deviation 


Figure  3.9:  Precision  of  range  estimate  is  shown.  Cramer-Rao 
lower  bounds  are  also  illustrated. 
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Figure  3.10:  (a)  Estimation  of  Gi  as  G2  varies  parametric  in  bias, 

(b)  Estimation  of  G2  as  G2  varies  parametric  in  bias. 
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Bias  Est.  (G^=2  p^=10ns  f^=10GHzt^=0.1ps) 
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Figure  3.11:  (a)  Bias  estimation  with  bias  py^  =  10  ns  parametric  in  bias, 

(b)  Bias  estimation  with  pulse  width  B  =  3  parametric  in  pulse  width. 


However,  the  bias  estimate  used  in  this  simulation,  Eq.  3.109,  omits  laser  pulse  data. 
Therefore,  the  sample  variance  of  the  bias  estimates  in  Fig.  3.11a  do  not  depend  on 

G2. 


In  Fig.  3.11b,  bias  estimation  is  simulated  for  different  values  of  G2  and  pulse 
width  (piu).  The  bias  is  fixed  at  H  =  3  in  all  simulations.  As  in  Fig.  3.11a,  the  CRLB 
is  a  function  of  more  than  just  the  bias.  In  this  case,  it  drops  as  pulse  width  decreases 
because  there  is  more  background-only  data  from  which  to  measure  the  bias.  And  as 
with  the  Fig.  3.11a  case,  the  CRLB  varies  with  G2.  The  simulations,  on  the  other 
hand,  do  not  use  the  laser  pulse  data  and  have  roughly  the  same  sample  variance  for 
all  levels  of  gain. 


3.3  Variable  Pulse  Width  LADAR  Signals 

The  single  parabolic  pulse  model  from  [4]  (discussed  in  Sec.  3.1)  and  the  ob¬ 
scured  target  model  from  Sec.  3.2  assume  that  the  pulse  width  in  the  received  signal  is 
a  known  constant.  This  assumption  is  valid  in  some  cases.  The  width  of  the  transmit¬ 
ted  pulse  in  a  LADAR  system  is  almost  always  known  to  a  high  degree  of  precision. 
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Therefore,  if  the  reflection  process  does  not  alter  the  pulse  width  by  a  signihcant 
amount,  then  the  received  pulse  width  is  the  same  as  the  transmitted  width. 

There  are  scenarios  where  the  reflection  process  does  not  make  a  significant 
change  to  the  pulse  width.  When  the  transmitted  laser  pulse  hits  a  surface  that  is 
perpendicular  to  the  LADAR  line  of  sight,  the  width  of  the  received  pulse  is  equal  to 
the  width  of  the  transmitted  pulse.  Depending  on  the  transmitted  pulse  width  and 
the  size  of  the  beam  after  propagation,  small  perturbations  from  perpendicular  in 
surface  orientation  may  not  change  the  pulse  length  by  a  significant  amount.  When 
the  beam  hits  a  surface  that  has  depth  in  range,  there  is  always  an  increase  in  the 
width  of  the  pulse  that  is  received.  But  if  the  target’s  range  depth  is  small  compared 
to  the  length  of  the  pulse,  then  the  increase  in  received  pulse  width  will  be  negligible. 

While  there  are  situations  where  the  received  laser  pulse  has  the  same  width 
as  the  transmitted  pulse,  there  are  also  many  real  world  scenarios  where  the  received 
LADAR  signal  is  wider  than  that  pulse  that  was  transmitted.  The  laser  pulse  can  hit 
a  flat  surface  that  is  not  perpendicular  to  the  line  of  sight.  In  that  case,  the  received 
pulse  is  stretched.  Objects  being  imaged  such  as  buildings,  vehicles,  and  vegetation 
usually  have  range  depth  that  will  increase  the  width  of  the  received  pulse.  In  general, 
the  amount  of  pulse  width  distortion  increases  with  range  since  the  width  of  a  focused 
laser  beam  increases  with  distance.  The  problem  is  also  exacerbated  as  laser  pulse 
width  decreases.  Short  pulses  are  desirable  for  good  range  precision,  but  are  more 
easily  distorted  by  targets  with  range  depth. 

In  this  section,  a  single  parabolic  pulse  model  is  used  to  simulate  the  received 
signal.  That  finite  duration  parabola,  plus  a  bias  level  used  to  simulate  noise  sources, 
models  the  mean  value  of  the  received  signal.  As  in  previous  sections,  the  received 
signal  is  a  sequence  of  statistically  independent  Poisson  random  variables.  But  in  this 
case,  the  pulse  width  is  treated  as  an  unknown  parameter  that  must  be  estimated. 
However,  it  is  shown  through  Cramer-Rao  lower  bound  analysis  that  when  the  pulse 
width  is  unknown,  that  the  range  precision  of  the  system  does  not  necessarily  decrease. 
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3.3.1  Variable  Pulse  Width  LADAR  Signal  Model.  The  signal  model  used 
for  the  variable  pulse  width  case  is  equal  to  the  parabolic  pulse  model  from  Eq.  3.2. 
The  model  is 


/(4)  =  G 


1  - 


4  -  2R/ c 

Pw 


rect  I  +  B. 


(3.112) 


Mathematically,  this  model  is  equivalent  to  the  one  presented  in  Sec.  3.1.  But  in  this 
case,  the  pulse  width  (p^)  is  an  unknown  that  must  be  estimated  from  the  data.  The 
unknown  pulse  width  will  increase  the  size  of  the  Fisher  information  matrix. 


3.3.2  Log-Likelihood  Function  and  its  Derivatives.  Because  the  model  in 
Eq.  3.112  has  an  unknown  parameter  that  is  not  present  in  previous  models,  the 
CRLBs  derived  for  those  models  can  not  be  used.  The  Fisher  information  matrix  for 
this  model  contains  elements  that  relate  to  the  unknown  pulse  width  p^,. 

The  log-likelihood  function  for  the  variable  pulse  width  LADAR  signal  model  is 

K  K  K 

l{R,G,B,p^)  =  ^ci(4)log[/(4)]  -  ^/(4)  -  ^  log  [(i(4)!] .  (3.113) 

k=l  k=l  k=l 

Derivatives  of  the  log-likelihood  function  are  used  to  calculate  the  CRLBs.  However, 
many  of  the  needed  derivatives  were  already  calculated  in  Sec.  3.1. 

The  first  derivatives  of  the  log-likelihood  function  with  respect  to  range,  gain, 
and  bias  have  already  been  computed.  The  derivative  with  respect  to  pulse  width  is 
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(3.114) 


The  only  second  derivative  diagonal  element  from  the  FIM  that  has  not  been  com¬ 
puted  previously  is 
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The  off-diagonal  elements  of  the  FIM  associated  with  pulse  width  are 
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The  derivatives  of  Eq.  3.112  with  respect  to  pulse  width  is  needed  to  calculate 
the  Fisher  information  matrix  elements  for  the  variable  pulse  width  model.  This 
derivative  is 
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The  elements  of  the  Fisher  information  matrix  are  the  negatives  of  the  expec¬ 
tations  of  the  second  derivatives  of  the  log-likelihood  function.  Many  of  these  have 
already  been  calculated  in  Sec.  3.1.  The  FIM  terms  that  are  needed  for  the  variable 
width  model  are 
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3.3.3  Fisher  Information  Matrix.  The  Fisher  information  matrix  for  the 
variable  pulse  width  model  from  Eq.  3.112  is 


(3.124) 
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First,  the  diagonal  elements  from  the  FIM  are  evaluated.  The  range,  gain,  and 
bias  were  already  derived  in  Sec.  3.1. 
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The  final  FIM  diagonal  element  is  the  one  associated  with  pulse  width.  Eqs.  3.120 
and  3.119  are  combined  to  obtain  the  following  integral  approximation. 
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where  the  change  of  variable  from  Eq.  3.26  was  used.  This  equation  can  be  evaluated 
using  the  integral  from  Eq.  3.40. 
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where  a  is  defined  in  Eq.  3.32.  This  analytic  approximation  is  compared  to  a  sum¬ 
mation  in  Fig.  3.12a.  Both  numbers  are  large  and  they  are  almost  identical. 

As  in  previous  cases,  odd  symmetry  causes  several  FIM  elements  to  be  ap¬ 
proximately  zero.  Inspection  of  the  hrst  derivative  equations  shows  that  the  range 
derivative  is  odd  about  the  pulse  center,  but  that  all  other  derivatives  are  even  about 
the  pulse  center.  When  the  terms  are  multiplied,  range  and  all  others  produce  odd 
functions  that  are  about  zero  after  integration.  Therefore, 
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The  odd  symmetry  of  the  function  associated  with  range  and  pulse  width  is  illustrated 
in  Fig.  3.12b.  The  functions  for  the  other  two  terms  were  present  in  previous  signal 
models  and  are  illustrated  in  Fig.  2.4.  This  range  and  pulse  width  function  was 
summed.  The  numeric  result  is  small  (0.06),  but  not  exactly  zero  since  the  sampling 
was  not  aligned  exactly  with  the  pulse  center. 


The  off-diagonal  FIM  element  relating  gain  and  bias  was  evaluated  earlier  (see 
Eq.  3.47). 
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The  FIM  element  relating  gain  and  pulse  width  must  be  calculated.  Eqs.  3.23,  3.119, 
and  3.122  can  be  combined  with  the  change  of  variable  from  Eq.  3.26  to  produce  the 


83 


-E[3\(R,G,B,pJ/a  Gd  pj  -E[3\(R,G,B,pJ/3  p^ 


.X  10 


The  3  function 


20 


40  60 

Time  (ns) 


(a) 


R=10 

G=5 

p^=20ns 

t^=100ns 

Sum:  637 
Approx.: 

992693.055 

537961811.5 

017 

\ 

V.  y 

80 


100 


X  10 


The  3  R3  p  function 


0.8- 

0.6 

?  0.4 
cn 

2,  0.2 

a.^ 

cq'  0 
O 

CC  -0.2 

S  -0-4 
-0.6- 
-0.8- 
-1 

0 


20 


/ 

R=10 

.....G-5 . 

/ 

B=2 

D  =20ns 

y 

t^=100ns 

Slim:  0.0 
Approx.: 

32009 

} 

40  60 

Time  (ns) 


80 


100 


(b) 


(c) 


(d) 


Figure  3.12:  (a)  The  pulse  width  second  derivative,  (b)  The  range-width  derivative, 

(c)  The  gain-width  derivative,  (d)  The  bias-width  derivative. 
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following  integral  approximation, 
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Using  Eqs.  3.29  and  3.40,  the  approximation  becomes 
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This  approximation  is  compared  to  a  numeric  integration  in  Fig.  3.12c.  The  results  are 
nearly  equal.  The  hnal  FIM  term  in  the  variable  pulse  width  model  is  a  combination 
of  Eqs.  3.123  and  3.119. 
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This  integral  is  proportional  to  the  one  that  was  evaluated  to  find  the  range  diagonal 
element  of  the  Fisher  information  matrix  (Eq.  3.25).  Therefore, 
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where  Eq.  3.88  was  used.  This  approximation  is  compared  to  a  numeric  integration 
in  Fig.  3.12d.  The  results  are  equal. 
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3.3.4  Cramer- Rao  Lower  Bounds. 
Fisher  information  matrix  can  be  written 


Using  the  results  from  Sec.  3.3.3,  the 
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This  is  a  block  diagonal  matrix.  It  can  be  inverted  using  Eq.  3.57.  The  Cramer-Rao 
lower  bound  for  range  is  the  reciprocal  of  the  upper-left  element  of  Eq.  3.136. 
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The  range  CRLB  is  identical  to  the  one  derived  for  the  original  parabolic  pulse  model 
(Eq.  3.61)  in  spite  of  the  fact  that  pulse  width  was  known  in  that  model  but  is 
unknown  in  this  case.  However,  range  estimation  is  more  difficult  when  the  width  is 
unknown. 

The  non-zero  3x3  sub-matrix  in  the  FIM  makes  the  CRLBs  for  gain,  bias, 
and  pulse  width  complicated.  Instead  of  presenting  complicated  expressions  for  those 
CRLBs,  they  are  listed  in  terms  of  Fisher  information  matrix  elements.  All  three 
CRLBs  calculations  use  the  determinant  of  the  sub-matrix  from  Eq.  3.136. 
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Using  this  result,  the  CRLBs  are 
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The  Fisher  information  matrix  elements  for  the  variable  pulse  width  model  from 
Eq.  3.112  are  summarized  in  Table  3.4.  The  CRLBs,  written  in  terms  of  the  FIM 
elements  are  shown  in  Table  3.5. 

3.3.5  Simulations.  Estimates  of  range,  gain,  and  bias  have  been  discussed 
and  tested  in  Secs.  3.1  and  3.2.  Though  the  range  CRLB  is  the  same  for  all  three  signal 
models,  the  actual  process  of  range  estimation  differs  depending  on  the  model.  In  this 
section,  estimation  of  range  and  pulse  width  are  simulated.  Parameter  estimation  for 
the  variable  pulse  width  model  from  Eq.  3.112  is  difficult  because  there  are  four 
unknown  parameters  that  must  be  estimated.  A  relatively  complicated  algorithm 
to  estimate  multiple  signal  parameters  is  derived  in  Sec.  5.2.2  and  tested  on  real 
LADAR  data.  In  this  section,  a  simpler  clairvoyant  pulse  width  estimation  algorithm 
is  developed  for  the  purpose  of  validating  the  width  CRLB.  The  process  of  range 
estimation  is  also  simulated  and  compared  to  the  range  CRLB. 
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Table  3.4:  Variable  Width  Model 

Fisher  Information  Matrix  Elements. 
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Table  3.5:  Variable  Width  Model 

Cramer-Rao  Lower  Bounds. 
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Figure  3.13:  Range  estimation  using  a  clairvoyant  (known 

received  pulse  width)  filter  and  a  filter  with  a  hxed  10  ns  hlter. 


The  range  estimation  simulations  from  Secs.  3.1  and  3.2  assumed  that  the  pulse 
width  in  the  received  signal  is  known.  In  practice,  the  width  may  be  unknown.  The 
simulation  in  Fig.  3.13  shows  the  precision  of  range  estimation  with  knowledge  of  the 
received  pulse  width  (clairvoyant)  and  without  knowledge  of  the  pulse  width. 

In  the  simulations  in  Fig.  3.13,  the  range  estimation  method  from  Sec.  3.1.4  was 
used.  Pulse  widths  between  10  and  20  ns  were  simulated.  In  the  clairvoyant  case,  the 
width  was  known  during  estimation.  In  the  other  case,  the  pulse  width  used  by  the 
hlter  was  hxed  at  10  ns  regardless  of  actual  pulse  width.  The  results  show  that  using 
the  wrong  pulse  width  in  the  hlter  will  decrease  the  precision  of  range  estimates. 

The  problem  of  pulse  width  estimation  is  now  addressed.  The  Gauss-Newton 
method  of  parameter  estimation  is  described  in  [23].  This  process  assumes  knowledge 
of  the  probability  distribution  of  the  data.  It  requires  an  initial  estimate  of  the 
unknown  parameter.  That  initial  estimate  is  updated  through  multiple  iterations 


until  it  converges  to  the  parameter  estimate.  The  iteration  is  a  function  of  the  data’s 
log-likelihood  function. 

The  unknown  signal  parameter  is  denoted  6.  The  initial  estimate  is  6q.  The 
Gauss-Newton  iterative  estimate  is 


9 


n 


902  90 


(3.142) 


The  Gauss-Newton  algorithm  must  be  initialized  with  a  coarse  estimate  of  the 
pulse  width.  A  clairvoyant  width  estimate,  that  is  an  estimate  that  assumes  knowledge 
of  the  gain,  range,  and  bias  of  the  signal,  can  be  derived  using  the  expectation  of  the 
sum  of  the  data.  The  expectation  of  the  sum  of  the  data  is 
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where  the  change  of  variable  u  =  tc  —  Rrt  was  used.  Thus,  the  following  equation  is 
an  unbiased  pulse  width  estimate. 
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This  becomes  the  initial  estimate  for  the  Gauss-Newton  iteration. 
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The  Gauss-Newton  iteration  is  derived  using  the  log-likelihood  function  asso¬ 
ciated  with  the  model  from  Eq.  3.112.  The  first  and  second  derivatives  of  the  log- 
likelihood  with  respect  to  pulse  width  were  already  computed.  They  are  shown  in 
Eqs.  3.114  and  3.115.  The  derivatives  of  the  mean  photocurrent  with  respect  to  pulse 
width  are 
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Therefore,  the  terms  from  the  Gauss-Newton  iteration  (Eq.  3.142)  are 
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When  evaluating  Eqs.  3.148  and  3.149,  the  current  pulse  width  estimate  is  used  to 
calculate  the  mean  photocurrent  I{tk)-  These  two  equations  are  used  to  compute  the 
Gauss-Newton  iteration  from  Eq.  3.142.  The  iteration  is  done  repeatedly  until  the 
pulse  width  estimate  converges  to  a  stable  value. 

The  Gauss-Newton  iteration  for  pulse  width  estimation  was  simulated  using 
Monte  Garlo  methods.  Many  realizations  of  Poisson-distributed  signals  were  gener¬ 
ated  for  different  gain,  bias,  and  pulse  width  values.  The  algorithm  was  initialized 
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Figure  3.14:  (a)  Pulse  width  estimation  as  width  and  gain  vary,  (b)  Pulse  width 

estimation  as  width  and  bias  vary. 

using  Eq.  3.145.  Twenty  iterations  of  the  algorithm  were  implemented  per  realization 
because  tests  showed  that  it  consistently  converged  after  that  number. 

The  results  of  the  Monte  Carlo  simulations  are  shown  in  Fig.  3.14.  In  Fig. 
3.14a,  the  bias  is  hxed  and  different  values  of  gain  are  simulated.  In  Fig.  3.14b,  the 
gain  is  constant  and  the  bias  is  varied.  Pulse  width  values  between  10  and  25  ns  are 
tested. 

Two  different  CRLBs  are  shown  in  Fig.  3.14.  The  solid  lines  describe  the  CRLB 
from  Eq.  3.141.  This  CRLB  does  not  assume  that  the  other  signals  parameters 
(range,  gain,  and  bias)  are  known.  The  simulations,  on  the  other  hand,  assumed  that 
those  parameters  were  known.  The  simulations  used  that  information  to  improve  the 
precision  of  the  pulse  width  estimate.  Because  the  clairvoyant  pulse  width  estimate 
knew  the  gain  and  bias,  that  estimate  was  more  precise  than  predicted  by  the  CRLB 
from  Eq.  3.141. 

The  second  CRLB  shown  in  Fig.  3.14  is  the  clairvoyant  pulse  width  estimate 
CRLB.  This  bound  assumes  knowledge  of  the  signal’s  gain  and  bias.  It  is  lower 
than  the  CRLB  from  Eq.  3.141.  That  CRLB  assumes  that  the  gain  and  bias  must 
be  estimated  from  the  data,  and  therefore,  will  not  be  as  precise  as  the  clairvoyant 
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estimate  that  knows  the  gain  and  bias.  The  clairvoyant  CRLB  is 


1 


This  CRLB  has  an  interesting  characteristic.  This  bonnd  increases  as  the  pnlse  width 
increases.  Eqnivalently,  the  CRLB  is  increasing  as  the  mean  pulse  energy  increases. 
In  most  cases,  an  increase  in  pulse  energy  implies  an  improvement  in  the  precision  of 
signal  parameter  estimates.  But  in  this  case,  if  the  pulse  energy  increase  is  due  to  an 
increase  in  pulse  width  while  the  gain  is  hxed,  the  pulse  energy  estimate  becomes  less 
precise. 

The  clairvoyant  pulse  width  CRLB  from  Eq.  3.150  is  always  less  than  or  equal 
to  the  pulse  width  CRLB  from  Eq.  3.141. 
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This  is  because  the  CRLB  for  any  parameter  is  always  greater  than  or  equal  to  the 
reciprocal  of  its  Fisher  information  matrix  element.  That  is. 
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See  Appendix  C  for  proof  of  this  property.  Thus  the  clairvoyant  CRLB  from  Eq. 


3.150  is  a  lower  bound  for  the  variance  of  unbiased  pulse  width  estimates.  However, 
it  is  not  as  tight  as  the  CRLB  from  Eq.  3.141. 

The  simulated  clairvoyant  estimates  shown  in  Fig.  3.14  of  pulse  width  are 
less  precise  than  the  clairvoyant  estimate  CRLB  says  they  could  be.  Gauss-Newton 
estimates  are  generally  not  efficient  [54].  That  is,  Gauss-Newton  estimates  generally 
have  a  variance  greater  than  the  Gramer-Rao  lower  bound.  However,  the  simulated 
clairvoyant  estimates  have  variance  lower  than  the  GRLB  from  Eq.  3.141.  That 
GRLB  assumes  the  estimator  does  not  know  the  gain  and  bias  and  is  higher  than  the 
clairvoyant  bound  from  Eq.  3.150. 
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3-4  Chapter  Summary 

In  this  chapter,  three  different  models  for  Poisson-distributed  LADAR  signals 
were  presented.  The  process  of  parameter  estimation  was  discussed  for  each  of  these 
models.  Cramer-Rao  lower  bounds  were  derived  for  all  parameters  present  in  these 
models  such  as  range  to  target,  gain,  and  bias. 

In  Sec.  3.1,  a  LADAR  signal  model  that  describes  a  laser  pulse  as  a  truncated 
parabolic  shape  was  discussed.  It  was  shown  that  a  relatively  compact  analytic  ex¬ 
pression  could  be  derived  for  that  model’s  range  CRLB.  In  Sec.  3.2,  a  model  with 
two  different  LADAR  signals  was  discussed.  This  model  is  used  to  describe  the  signal 
observed  when  a  target  is  partially  obscured  by  forest  canopy  or  camouflage  netting. 
It  was  shown  that  the  presence  of  the  second  signal  did  not  decrease  the  range  CRLB 
to  the  target.  Finally,  in  Sec.  3.3,  the  effect  of  unknown  pulse  width  was  examined.  It 
was  shown  that  lack  of  knowledge  of  the  pulse  width  did  not  affect  the  range  CRLB. 
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IV.  Range  Precision  of  LADAR  Systems 

The  topic  of  this  chapter  is  the  precision  of  range  estimates  made  by  LADAR  systems. 
In  the  previous  chapter,  it  was  assumed  that  the  noise  in  a  LADAR  signal  follows 
the  Poisson  distribution  and  that  the  shape  of  a  received  laser  pulse  was  described 
by  a  truncated  parabola.  That  model  is  discussed  in  this  chapter,  but  other  methods 
of  signal  modeling  are  also  used.  In  some  cases,  the  noise  in  the  LADAR  system 
is  assumed  to  have  a  Gaussian  distribution.  The  shape  of  the  laser  pulse  can  be 
modeled  using  functions  other  than  the  truncated  parabola.  In  this  chapter,  the  use 
of  alternative  pulse  shape  models  is  discussed. 

In  Sec.  4.1,  the  Poisson-distributed  parabolic  pulse  model  from  Chapters  II 
and  III  is  employed  to  simulate  shot-noise  limited  LADAR  signals.  The  precision  of 
range  estimates  made  with  shot-noise  limited  systems  is  discussed.  The  effect  of  pulse 
averaging  is  explored  and  the  results  are  compared  to  systems  with  signals  that  have 
white  Gaussian,  rather  than  Poisson  noise. 

In  Sec.  4.2,  the  effect  of  target  surface  orientation  on  range  precision  is  quan- 
tihed.  It  is  shown  that  for  typical  LADAR  systems  this  effect  is  only  signihcant  for 
targets  that  are  many  kilometers  away.  In  order  to  simplify  the  analytics  involved  in 
this  problem,  the  laser  pulse  is  assumed  to  have  a  Gaussian  shape  and  the  noise  is 
modeled  as  a  white  Gaussian  process.  Cramer-Rao  lower  bounds  are  calculated  for 
range  precision  given  a  tilted  target  surface. 

4.1  Target  Range  Estimation 

A  Cramer-Rao  lower  bound  for  range  estimate  variance  was  derived  in  Chapter 
III.  Methods  of  range  estimation  were  also  developed  and  tested  in  that  chapter. 
It  was  noted  that  for  shot-noise  limited  (Poisson  distributed)  data,  the  maximum 
likelihood  estimate  could  not  be  solved  analytically.  However,  it  was  also  shown 
that  a  combination  of  matched  hltering  and  peak  fitting  produced  reasonably  precise 
range  estimates  in  the  sense  that  they  were  close  to  the  CRLB.  However,  it  was  also 
demonstrated  when  the  matched  filter  was  replaced  by  a  filter  equal  to  the  square 
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root  of  the  laser  pulse  shape,  the  range  precision  improved  in  most  cases.  Range 
estimates  produced  by  the  square-root  hlter  nearly  achieved  the  CRLB  when  the 
signal  amplitude  was  large  compared  to  the  background  level  (see  Fig.  3.3). 

In  this  section,  the  Poisson  noise  range  estimation  CRLB  and  simulation  results 
are  compared  to  range  estimation  for  a  signal  mixed  with  additive  white  Gaussian 
noise  (AWGN).  It  is  shown  that  the  Gaussian-distributed  signal  produces  range  esti¬ 
mates  that  are  more  precise  than  can  be  obtained  from  a  Poisson-distributed  signal 
with  equal  background  noise  variance  and  mean  signal  energy.  This  is  proved  through 
simulation  and  by  comparing  the  Gramer-Rao  lower  bounds  for  these  signals. 

4.1.1  Comparison  to  Gaussian-Noise  Signals.  Noise  in  radar  and  commu¬ 
nication  systems  is  often  assumed  to  have  a  Gaussian  distribution.  Gaussian  noise 
has  also  been  used  to  model  LADAR  signals  [9,31].  Gaussian  noise  models  are  com¬ 
mon  in  the  communication  and  signal  processing  literature  and  are  well-understood. 
For  purposes  of  comparison,  the  Poisson  noise  range  GRLB  derived  in  Ghapter  III 
is  compared  to  a  Gaussian  noise  range  GRLB  that  is  derived  in  this  section.  The 
comparison  uses  the  parabolic  pulse  signal  model  from  Eq.  3.112  for  both  Poisson 
and  Gaussian  noise. 

The  signal  model  is 


(4.1) 


It  is  assumed  that  Gaussian-distributed  samples  dg{tk)  have  been  collected.  The  mean 
value  of  these  samples  is  the  mean  photocurrent  I{tk).  The  variance  of  the  samples 
is  B  since  that  is  the  variance  of  the  background  in  the  Poisson-distributed  case. 


113(4)  ~  AT  [/(4),  5]. 


(4.2) 


96 


As  in  the  Poisson  case,  it  is  assumed  that  all  samples  are  statistically  independent  of 
each  other.  Note  that  in  this  case,  it  is  possible  for  the  signal  to  have  negative  values. 

The  range  estimate  CRLB  for  a  signal  mixed  with  AWGN  was  derived  in  Fun¬ 
damentals  of  Statistical  Signal  Processing:  Estimation  Theory  by  Kay  [23] . 
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where  Rg  is  any  unbiased  range  estimate.  This  Gaussian-noise  model  GRLB  derived 
in  [23]  was  calculated  assuming  that,  in  the  absence  of  signal,  the  noise  is  zero-mean. 
It  is  assumed  here  that  there  is  a  non-zero  background  mean  of  B.  However,  the 
Gaussian-noise  range  GRLB  is  unaffected  by  the  presence  of  a  bias  (the  background 
non-zero  mean)  in  the  data.  For  more  detail  on  why  the  Gaussian-noise  range  GRLB 
is  not  dependent  on  bias,  see  [16]. 

The  integral  in  Eq.  4.3  is 
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where  the  change  of  variable  v  =  t  —  2R/c  was  used.  Therefore,  the  Gaussian-noise 
range  GRLB  is 
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For  the  purpose  of  comparison,  a  bound  on  the  Poisson-noise  range  GRLB  (Eq. 
3.137)  is  desired.  That  GRLB  is  (see  Eqs.  3.137  and  3.32) 
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The  bound  will  be  used  to  compare  this  GRLB  to  the  Gaussian-noise  GRLB  from  Eq. 

4.5. 
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Inverse  hyperbolic  tangent  can  be  written  as  an  infinite  series  [1]. 
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where  |x|  <  1.  If  atanh(x)  is  divided  by  x,  then  one  is  subtracted,  the  geometric 
series  can  be  used  to  find  an  upper  bound.  Consider  the  following  series  expansion. 
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Using  the  fact  that 
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Eq.  4.8  can  be  written 
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Therefore,  the  term  from  the  denominator  of  the  Poisson-noise  CRLB  (Eq.  4.6) 
is  bounded  by 

v/^atanh  -  1  <  ^.  (4.11) 

This  implies  that  the  Poisson-distribution  CRLB  is  bounded  by 
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This  lower  bound  on  the  Poisson-distribution  CRLB  is  equal  to  the  Gaussian-noise 
CRLB  from  Eq.  4.5.  If  the  CRLB  can  be  achieved,  then  the  variances  of  the  estimators 
are  related  by 
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Parabolic  Pulse  Range  Est.  (B=10  K=100) 


Parabolic  Pulse  Range  Est.  (P^^=‘I0  K=100) 


Figure  4.1:  (a)  Background  is  i?  =  10  and  pulse  width  varies. 

(b)  Pulse  width  is  pw  =  10  ns  and  background  varies. 

The  CRLB  is  sometimes  used  to  design  LADAR  systems.  Eq.  4.12  shows  that 
the  Gaussian-noise  CRLB  is  always  less  than  the  Poisson  signal  CRLB.  The  Poisson 
CRLB  is  therefore  a  more  realistic  representation  of  the  range  precision  of  a  shot  noise 
limited  LADAR.  Because  the  Gaussian  CRLB  is  always  less  than  the  Poisson  CRLB, 
Eq.  4.5  is  a  lower  bound  for  shot  noise  limited  range  precision.  However,  the  Poisson 
CRLB  (Eq.  4.6)  is  a  tighter  bound  on  the  range  precision  limit  of  a  shot  noise  limited 
LADAR  system.  That  is,  real  unbiased  range  estimates  will  have  variances  closer  to 
the  Poisson  CRLB  than  to  the  Gaussian  CRLB. 

4- 1-2  Simulations.  Monte  Carlo  simulations  of  range  estimation  were  com¬ 
pleted  and  compared  to  the  Poisson  and  Gaussian  CRLBs.  The  results  are  shown  in 
Fig.  4.1.  In  the  simulations,  a  =  10  ns  parabolic  pulse  was  sampled  at  frequency 
/s  =  1  GHz.  One  hundred  samples  were  simulated  {K  =  100),  so  the  sampling  dura¬ 
tion  was  td  =  100  ns  where  td  =  K/ fg.  The  bias  level  in  the  simulation  was  B  =  10. 
This  value  {B)  is  also  the  variance  of  the  background  samples  in  the  Poisson  and 
Gaussian  distribution  cases. 
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The  pulse  amplitude  (G)  was  varied  from  10“^  to  10^.  For  each  pulse  amplitude, 
10,000  Poisson  and  10,000  Gaussian  signal  realizations  were  generated.  The  mean 
value  of  these  signals  was  the  parabolic  pulse  from  Eq.  4.1.  Each  of  the  signal 
realizations  was  match  filtered.  The  simulated  range  estimates  were  calculated  using 
three-point  parabolic  fits  to  the  peaks  of  the  match  filter  outputs.  The  parabolic  fit 
equations  used  are  shown  in  [20]. 

In  Fig.  4.1,  the  CRLBs  from  Eqs.  4.5  and  4.6  are  compared  to  the  sample 
variances  produced  by  the  Monte  Carlo  simulations.  At  pulse  amplitudes  greater 
than  the  bias  level  (G  >  B),  the  bounds  are  tight  for  the  Poisson  and  Gaussian 
distribution  cases.  As  the  pulse  amplitude  drops  below  the  bias  level,  the  simulation 
sample  variances  begin  to  deviate  from  the  CRLBs.  At  very  low  signal  amplitudes 
(G  <<  B),  the  range  estimation  algorithm  is  reporting  a  random  variable  that  is 
uniformly  distributed  over  the  sampling  interval.  In  uniform  distribution  case,  the 
range  variance  is  Eq.  3.79. 

At  very  low  signal  amplitudes,  the  sample  variance  is  less  than  the  CRLB. 
Cramer-Rao  lower  bounds  are  only  valid  for  unbiased  estimators.  In  the  very  low 
amplitude  case,  the  range  estimate  is  biased.  Its  mean  value  is  near  the  center  of 
the  sampling  interval  regardless  of  the  laser  pulse’s  true  position.  It  is  this  bias  that 
allows  the  estimator’s  variance  to  surpass  the  bound.  But  in  this  case,  the  estimator’s 
mean  squared  error  (variance  plus  bias  squared)  can  be  very  high.  The  fact  that 
sample  variance  is  below  the  CRLB  should  obviously  not  be  mistaken  for  efficient 
performance. 

The  bias  in  low  amplitude  range  estimates  is  illustrated  in  Fig.  4.2.  The  bias  is 
shown  using  Monte  Carlo  simulations  where  one-hundred  sample  LADAR  signals  are 
simulated.  A  parabolic  pulse  is  centered  at  sample  70  and  the  bias  level  is  B  =  10.  In 
Fig.  4.2a,  the  gain  is  one  (G  =  1).  The  pulse  is  visible,  but  as  the  gain  is  so  low  that 
the  distribution  is  nearly  uniform  over  the  sampling  interval.  The  range  estimates 
are  biased  toward  the  center  of  the  sampling  interval.  The  average  range  estimate  is 
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Histogram  of  Range  Estimates 


Estimated  Range 

(a) 


Histogram  of  Range  Estimates 


Estimated  Range 

(b) 


Figure  4.2:  (a)  Histogram  of  range  estimates  with  gain  G  =  1. 

(b)  Range  estimates  histogram  with  gain  G  =  3. 


56.9.  In  Fig.  4.2b,  the  range  gain  has  increased  to  three  (G  =  3).  Outliers  are  visible 
throughout  the  sampling  interval.  However  outliers  are  rare  and  do  little  to  bias  the 
range  estimate.  The  average  estimate  of  67.7  is  close  to  the  true  value. 


4-.  1.3  Multiple  Pulses  and  Range  Precision.  Range  estimates  from  multiple 
LADAR  signals  can  be  combined  to  improve  precision.  Multiple  range  estimates  could 
be  averaged,  or  signals  could  be  accumulated  before  one  range  estimate  is  computed. 
Regardless  of  the  method  of  processing,  (e.g.  averaging  of  estimates  or  accumulation 
then  estimation),  the  Cramer- Rao  lower  bound  on  range  precision  is  valid  for  any 
unbiased  range  estimate. 

Given  N  statistically  independent  signals,  the  variance  CRLB  for  any  parameter 
estimate  decreases  by  a  factor  of  N  [5].  In  that  sense,  having  multiple  pulses  helps 
range  precision.  However,  this  assumes  that  the  laser  has  transmitted  A^-times  more 
energy  than  was  transmitted  in  a  single  pulse. 

If  the  amount  of  laser  energy  is  hxed,  but  can  be  divided  equally  into  N  pulses, 
then  the  CRLB  varies  with  N.  In  a  photon  counting  system,  mean  pulse  energy  is 
proportional  to  the  mean  photocount  gain  G  (see  Eqs.  2.26  and  2.27).  The  range 
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CRLB,  given  that  a  total  (accumulated)  amplitude  of  G  is  divided  equally  into  N 
pulses  is 
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where  R  denotes  any  unbiased  range  estimate  obtained  from  N  Poisson-distributed 
pulses,  each  with  amplitude  G/N. 


In  order  to  compare  this  CRLB  to  the  single  pulse  CRLB  from  Eq.  4.6,  the 
inhnite  series  for  inverse  hyperbolic  tangent  (see  Eq.  4.10)  is  used. 
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Gain  G  and  background  B  are  always  positive.  Therefore, 
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This  inequality  implies  that 
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Because  the  number  of  pulses  is  a  positive  integer. 
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with  equality  if  and  only  if  =  1.  This  equation  means  that  the  CRLBs  from  Eq. 
4.6  and  4.14  are  related  by 
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If  range  estimators  can  be  fonnd  that  achieve  the  Cramer-Rao  lower  bonnd,  then  the 
CRLB  for  a  single  Poisson-distribnted  laser  pnlse  range  measnrement  R  and  a  range 
measnrement  made  from  a  set  of  N  Poisson-distribnted  pulses  R  are  related  by 


Var 


>  Var 


(4.20) 


with  equality  if  and  only  if  iV  =  1. 

Eq.  4.20  proves  that  for  a  hxed  total  laser  energy  and  Poisson-distribnted 
signals,  the  range  precision  is  maximized  when  all  energy  is  transmitted  in  a  single 
pulse.  It  is  interesting  to  compare  this  fact  to  Poisson-distribnted  signal  detection 
results  from  Gatt  and  Henderson  [9].  In  that  paper,  it  was  shown  that  for  a  fixed 
total  energy,  the  detectability  of  Poisson  signals  is  maximized  by  transmitting  all 
energy  in  a  single  pulse.  Thus,  there  is  a  relationship  between  the  range  precision  and 
the  detectability  of  a  Poisson-distribnted  signal. 

In  [35],  the  CRLB  for  position  estimates  of  a  Gaussian  pulse  in  Poisson  noise 
was  calculated.  In  that  paper,  it  was  assumed  that  there  was  no  background  noise 
{B  =  0).  The  result  in  Eq.  4.20  stands  in  contrast  to  the  zero  background  noise  result 
(Eq.  18  of  [35]).  That  CRLB  is  inversely  proportional  to  amplitude.  Thus,  the  range 
precision  is  independent  of  the  number  of  averages  in  the  zero-background  noise  case 
when  the  total  laser  energy  is  hxed. 

The  CRLB  for  a  shot-noise  limited  system  which  divides  a  hxed  total  energy 
into  N  pulses  (Eq.  4.14)  is  illustrated  in  Fig.  4.3.  For  the  case  considered  in  the 
illustration  (p^  =  10  ns,  /s  =  1  GHz,  B  =  10,  and  G  varying  between  10“^  and  10^) 
the  range  precision  decreases  rapidly  as  laser  energy  is  divided  into  multiple  pulses. 
For  example,  when  an  accumulated  gain  of  ten  {G  x  N  =  10)  is  divided  into  three 
pulses,  the  range  standard  deviation  CRLB  is  about  twice  as  high  as  if  all  energy  was 
transmitted  in  a  single  pulse  over  most  accumulated  gain  values  shown  in  the  hgure. 
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Figure  4.3:  (a)  Effect  of  pulse  averaging  on  CRLB  with  B  =  1. 

(b)  Effect  of  pulse  averaging  with  B  =  10. 

4-2  Surface  Orientation  and  Range  Estimation 

LADAR  systems  are  capable  of  transmitting  laser  pulses  with  very  short  du¬ 
rations.  Typical  systems  have  pulses  that  are  only  a  few  nanoseconds  in  duration. 
Some  future  systems  will  employ  pulses  with  even  shorter  durations.  These  short 
pulses  provide  excellent  range  resolution  and  precision.  However,  when  a  target  has 
depth  over  the  area  illuminated  by  the  beam,  the  reflection  process  causes  the  received 
laser  signal  to  be  wider  than  the  transmitted  signal.  The  widened  received  signal  has 
less  range  precision  than  a  signal  that  is  not  elongated  due  to  target  depth.  Even 
when  a  flat  surface  is  targeted,  that  surface  has  depth  if  it  is  tilted  with  respect  to 
the  LADAR  line-of-sight. 

In  this  section,  the  loss  in  range  precision  due  to  the  tilt  of  a  flat  surface  is 
quantihed.  Unlike  in  previous  sections  of  this  dissertation,  the  noise  in  the  LADAR 
signal  is  assumed  to  be  Gaussian.  In  practice,  the  Gaussian  noise  model  does  not 
completely  capture  the  fluctuations  in  signal  level  due  to  laser  speckle.  However  it  is  an 
accurate  approximation  for  any  system  with  a  signihcant  amount  of  electronic  noise  or 
for  a  LADAR  that  makes  use  of  diversity  to  reduce  speckle  fluctuations.  The  temporal 
shape  of  the  pulse  is  also  modeled  as  a  Gaussian  function.  These  assumptions  simplify 
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the  analytic  calculations  need  to  analyze  this  problem  and  produce  concise  formulas 
describing  the  effect  of  tilt  on  the  range  CRLB. 

The  problem  of  target  depth  distorting  the  shape  of  received  laser  pulses  has 
been  addressed  by  other  authors  [15,45,46].  The  novel  material  contained  in  this 
section  is  the  analytic  expression  derived  for  a  received  laser  pulse  and  the  Cramer- 
Rao  lower  bound  on  range  precision  when  the  target  is  tilted. 

4.2.1  Scenario  and  Assumptions.  A  TEMOO  Gaussian  beam  is  assumed  to 
have  been  transmitted.  A  short  pulse  is  used  to  achieve  good  range  precision.  Laser 
pulses  are  described  by  their  duration  and  their  size.  The  size  is  the  beam’s  width 
transverse  to  the  system’s  LOS.  Propagated  laser  beam  size  is  a  function  of  beam 
size  at  the  transmit  aperture,  laser  wavelength,  the  focus  distance,  and  distance  of 
propagation.  Beam  size  is  dependent  on  the  LADAR  transmitter  aperture.  If  the 
width  of  the  transmit  aperture  is  hxed,  then  as  beam  size  increases  more  laser  power 
is  truncated  by  the  opaque  portions  of  the  system.  If  the  beam  size  is  greater  than 
the  aperture  radius,  then  a  signihcant  amount  of  laser  energy  is  truncated.  Thus,  the 
beam  size  is  almost  always  smaller  than  the  aperture  radius. 

At  ranges  of  less  than  a  few  kilometers,  the  effect  of  target  tilt  on  range  preci¬ 
sion  is  negligible.  However,  LADAR  systems  are  sometimes  operated  at  long  range. 
Systems  are  used  to  do  three-dimensional  imaging  from  the  air.  Airborne  LADAR 
systems  are  used  to  image  terrain  and  structures  [43].  The  three-dimensional  data 
can  be  used  for  mapping  and  for  automatic  target  recognition  (ATR)  [53].  Airborne 
systems  have  a  wide  held  of  regard  and  can  operate  at  very  long  ranges.  As  range 
increases,  the  beam  size  gets  large  and  the  effect  of  target  tilt  on  range  precision 
becomes  noticeable. 

Because  of  size,  weight,  and  power  restrictions,  transmit  apertures  in  airborne 
LADAR  systems  are  usually  only  several  inches  wide  at  most.  These  relatively  small 
apertures  keep  the  beam  size  small.  Smaller  transmitted  beam  size  means  larger  beam 
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size  in  the  far  field.  Increases  in  beam  size  at  the  target  exacerbate  the  problem  of 
range  precision  loss  due  to  target  tilt. 

In  military  imaging  scenarios  the  range  to  the  target  is  often  much  greater  than 
the  platform’s  altitude.  In  this  situation,  even  flat  ground  has  a  large  tilt  angle  with 
respect  to  the  LOS.  The  loss  in  range  precision  can  be  signihcant  in  this  case. 

It  is  assumed  that  the  LADAR  is  either  monostatic  or  that  the  transmitter  and 
receiver  apertures  are  located  close  together.  Monostatic  systems  use  transmit-receive 
(TR)  switches  to  prevent  the  detector  from  being  damaged  by  the  relatively  intense 
transmitted  laser  pulses.  A  TR  switch  enables  the  system  to  use  a  single  aperture 
for  the  transmission  and  reception  of  laser  energy.  If  the  system  is  bistatic,  then 
the  transmit  and  receive  apertures  are  almost  always  located  so  close  together  (they 
commonly  share  a  gimballed  turret  in  airborne  systems)  that  they  can  be  treated  as 
monostatic  for  the  purpose  of  this  analysis. 

4-2.2  Scenario  Geometry.  The  scenario  considered  in  this  section  is  illus¬ 
trated  in  Fig.  4.4.  A  short  duration  laser  pulse  is  transmitted  along  the  LADAR’s 
line-of-sight.  The  temporal  prohle  is  assumed  to  be  a  Gaussian  shape.  The  beam  has 
an  irradiance  prohle  that  is  dehned  by  the  amount  of  power  transmitted  over  time 
and  the  beam’s  size  after  propagating  to  the  target. 

The  beam  hits  a  surface  that  may  be  tilted  with  respect  to  the  LADAR’s  LOS. 
The  surface  is  assumed  to  be  hat  over  the  width  of  the  irradiance  prohle.  This 
assumption  produces  relatively  simple  analytic  solutions  for  the  received  signal.  If 
the  surface  is  not  hat,  then  shape  of  the  rehected  pulse  could  be  complicated  and 
would  likely  have  to  be  calculated  numerically.  This  sort  of  analysis  has  been  done 
in  [45,46]. 

The  rehected  laser  pulse  is  shown  in  Fig.  4.5.  The  target  tilt,  which  is  pa¬ 
rameterized  by  angle  (j)  (measured  in  radians),  has  elongated  the  rehected  pulse  as 
compared  to  the  incident  pulse  shown  in  Fig.  4.4.  Precision  of  a  range  measurement 
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Figure  4.4:  Illustration  of  a  laser  pulse,  its  irradiance  pattern, 
and  a  tilted  surface. 


is  dependent  on  the  width  of  the  received,  not  the  transmitted  pulse.  Thus,  the  target 
tilt  with  respect  to  the  LADAR  LOS  causes  a  drop  in  range  precision. 


4.2.3  Propagation  and  Reflection.  The  held  at  the  transmit  aperture  is 
assumed  to  be  Gaussian  with  beam  size  Wq.  The  beam  is  focused  at  Fq.  The  beam 
size  is  the  radius  where  the  held  magnitude  is  below  the  peak  magnitude  by  a  factor 
of  e~^.  This  Gaussian  held  is  given  by  Eq.  2.11  and  is  shown  again  here. 
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Figure  4.5:  Reflected  laser  pulse  from  a  surface  tilted  at  angle 
(j)  with  respect  to  the  LADAR  line-of- sight. 


The  function  Pt{t)  is  the  instantaneous  transmitted  power  in  the  laser  beam.  The 
Fresnel  propagation  of  the  field  from  Eq.  4.21  is 
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(4.23) 


The  variable  z  is  the  distance  of  propagation  along  the  LADAR  LOS.  The  propagated 
irradiance  is 

Io{r,z)  =  \Uo{r,z)\^  =  exp  (-2^)  (4-24) 

where  the  beam  size  at  range  z  is 


W,  =  Wo 


(4.25) 


Eqs.  4.22  -  4.25  were  also  shown  in  Sec.  2.1.3,  though  the  instantaneous  power  was 
not  shown  as  a  function  of  time  in  that  section. 
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In  a  typical  scenario,  a  LADAR  receives  only  a  tiny  fraction  of  the  laser  energy 
that  it  transmits.  If  the  laser  beam  reflects  off  of  a  flat  surface  that  is  oriented 
perpendicular  to  the  line-of-sight,  then  the  signal  received  by  the  LADAR,  which  is 
denoted  Pr{t),  has  the  same  shape  as  the  transmitted  pulse.  Thus,  the  received  signal 
is  proportional  to  the  transmitted  signal. 


Pr{t)  oc  Pt{t). 


(4.26) 


Note  that  the  delay  due  to  the  round-trip  propagation  time  has  not  been  introduced 
into  Pr{t).  That  delay  is  added  later  in  this  section  when  the  effect  of  target  tilt  is 
analyzed. 

Eq.  4.26  assumes  that  the  propagated  beam  size  is  small  compared  to  the 
distance  of  propagation.  That  is,  «  z.  If  a  wide  beam  is  propagated  a  short 
distance,  then  the  received  signal  Pr{t)  could  have  a  shape  different  than  the  shape  of 
the  transmitted  signal  Pt{t)-  This  happens  when  the  LOS  distance  between  different 
parts  of  the  illuminated  target  varies  by  an  amount  that  is  signihcant  compared  to 
the  laser  pulse  duration. 

It  is  assumed  that  the  LADAR  receiver  aperture  is  located  at  the  origin  of  a 
Cartesian  coordinate  system.  The  transmitted  laser  beam  propagates  in  the  positive 
2;-direction.  Referring  to  Figs.  4.4  and  4.5,  the  distance  of  propagation  from  the  origin 
to  target  at  coordinate  {x,  y)  is  approximately 


z  ^  R  +  tan(0)?/ 


(4.27) 


where  R  is  the  range  from  the  receiver  to  the  target  along  the  axis  of  propagation. 
The  preceding  approximation  assumes  that  the  vertical  coordinate  is  small  compared 
to  the  target  range  {y  «  R).  In  the  figures,  the  y-a.xis  has  been  defined  as  the 
vertical.  However,  the  orientation  is  arbitrary  because  of  the  radial  symmetry  of  the 
beam’s  irradiance  pattern. 
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A  tilt  in  the  orientation  of  the  target  surface  (0  7^  0),  will  distort  the  received 
laser  pulse  shape  as  compared  to  the  transmitted  shape.  The  signal  received  from  a 
small  rectangle  with  widths  /S.x  and  A?/  located  at  (x,  y)  is 


AxA2//o(x,|/,r)Pr 

A  A  2  /  +  ?/2\  /  2R\ 


(4.28) 
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where  c  is  the  speed  of  light.  The  radial  coordinate  r  in  the  propagated  irradiance 
pattern  (Eq.  4.24)  has  been  replaced  with  the  Cartesian  coordinate  {x,y)  where 
=  x^  +  2/^.  The  distance  of  propagation  equation  (Eq.  4.27)  was  also  used  to  arrive 
at  this  result. 

The  total  received  signal  from  a  tilted  surface  Ptiit{t)  can  be  derived  by  summing 
over  all  of  the  signals  from  the  small  rectangles  (Eq.  4.28). 
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This  result  can  be  approximated  by  an  integral. 


CO  CO 


2  (  x^  +  y^' 

exp  —2 — —  1  Pr 


TlW^ 


1E2 


^  2R  2tan(0) 
c  c  ^ 


dydx.  (4.30) 


x=— CO  y=—oo 


Integrating  over  the  x-axis  yields 
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If  the  surface  tilt  is  zero  (0  =  0),  then  the  this  equation  is  just  Prit),  delayed 
by  the  round-trip  propagation  time.  If  the  surface  tilt  is  non-zero  (0  7^  0),  then  the 
following  two  changes  of  variable  are  used  to  evaluate  Eq.  4.31. 
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Applying  these  changes  of  variable  to  Eq.  4.31  produces 
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where  *  denotes  the  convolution  operation. 


(4.32) 

(4.33) 


v)  dv 


c 

2  tan(0) 


(4.34) 


4.2.4  Gaussian  Pulses.  The  equation  describing  the  shape  of  a  Gaussian 
pulse,  which  is  denoted  Pgit),  is 

Pg{t)  =exp(^-^t^^  (4.35) 

where  w  is  the  pulse  width.  The  full-width  half-maximum  (FWHM)  is  a  commonly 
used  metric  to  describe  the  width  or  duration  of  a  pulse.  The  full- width  half  maximum 
for  the  Gaussian  pulse,  which  is  denoted  Tg,  is 

rg  =  n;[81og(2)]'/^  (4.36) 
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The  variable  E  is  the  laser  pulse  energy  that  is  intercepted  by  the  LADAR’s 
receiver  aperture  and  focused  onto  the  detector.  The  energy  is 


OO 


(4.37) 


This  equation  describes  the  amount  of  optical  energy  in  the  received  signal.  The  signal 
Pr{t)  is  not  the  usual  voltage  or  current  signal  common  in  the  electrical  engineering 
literature.  Thus,  the  signal  is  not  squared  to  calculate  its  energy  as  is  normally  done 
when  characterizing  signals. 

The  equation  for  a  received  Gaussian  pulse  with  energy  E  is 


(4.38) 


This  is  Eq.  4.35  after  scaling  so  that  its  energy  is  E. 


The  received  signal  from  a  flat  tilted  surface,  which  is  calculated  by  substituting 
Eq.  4.38  into  Eq.  4.34,  is 
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The  convolution  of  two  Gaussian  waveforms  with  widths  wi  and  W2  and  unity  peak 
height  is 


exp 


*  exp 
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Thus,  the  convolution  of  these  two  Gaussian  functions  is  another  Gaussian  function 
with  width  Using  Eq.  4.40,  Eq.  4.39  becomes 


(4.41) 
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The  variable  Wt  is  introduced  to  describe  the  width  of  a  received  Gaussian  pulse  that 
has  reflected  off  a  tilted  surface.  It  is 


2  2  I 

W.  =  W  + 


(4.42) 


Using  Wt,  Eq.  4.41  becomes 
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Examples  of  received  Gaussian  pulses  after  reflection  off  a  tilted  surface  are 
shown  in  Fig.  4.6.  In  this  figure,  a  2.5  cm  beam  propagates  10  km.  The  wavelength 
is  1.6  /im.  The  beam  hits  targets  tilted  at  various  angles  with  respect  to  the  LADAR 
LOS. 

In  Fig.  4.6a,  a  =  1  ns  FWHM  pulse  is  transmitted.  For  a  15  degree  tilt  angle, 
the  received  FHWM  is  only  1.09  ns.  But  with  a  45  degree  tilt,  the  pulse’s  FWHM 
nearly  doubles.  When  the  tilt  reaches  75  degrees,  the  received  FWHM  is  over  six 
times  greater  than  the  transmitted  FHWM.  In  Fig.  4.6b,  the  transmitted  FWHM 
is  2  ns.  This  wider  pulse  is  affected  less  by  surface  orientation  than  the  1  ns  pulse. 
However,  target  tilt  still  affects  the  received  signal.  For  a  45  degree  tilt,  the  received 
FWHM  increases  about  28%.  For  the  75  degree  tilt,  the  FWHM  is  over  three  times 
wider  than  the  transmitted  FWHM. 
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Pulses  Reflected  off  Tilted  Surface  Pulses  Reflected  off  Tilted  Surface 


Figure  4.6:  (a)  Received  signals  for  1  ns  FWHM  pulse, 

(b)  Received  signals  for  2  ns  FWHM  pulse. 


4-2.5  Parabolic  Pulses. 
energy  E  is 


-  ip; 


The  equation  for  a  received  parabolic  pulse  with 
3E 
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This  is  Eq.  3.74,  scaled  so  that  its  energy  is  E.  Also,  Eq.  3.74  is  being  evaluated 
at  range  zero  {Rrt  =  0)  in  this  case.  For  the  parabolic  pulse,  the  FWHM,  which  is 
denoted  Tp,  is 

Tp  =  V2p^.  (4.45) 


The  received  signal  from  a  flat  tilted  surface,  which  is  calculated  by  substituting 
Eq.  4.44  into  Eq.  4.34,  is 
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This  equation  is  proportional  to  the  convolution  of  a  Gaussian  (Eq.  4.35)  and  a 
parabolic  (Eq.  3.74)  pulse.  The  constants  in  Eq.  4.46  are  ignored  for  the  moment 
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and  the  convolution  is  denoted  y{t). 
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where 
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After  extensive  algebra,  it  can  be  shown  that  this  convolution  is 
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Eq.  4.49  is  illustrated  and  validated  in  Fig.  4.7.  The  Gaussian  pulse  width 
is  hxed  at  ta  =  2  ns  and  the  parabolic  pulse  width  varies  between  0.5  and  4  ns. 
The  pulses  are  shown  in  Fig.  4.7a  and  the  convolution  in  Fig.  4.7b.  The  analytic 
expression  for  the  convolution  is  shown,  along  with  a  numeric  calculation  in  order  to 
validate  the  complicated  expression  from  Fig.  4.7. 

Eq.  4.49  can  be  used  in  Eq.  4.46  to  calculate  the  received  signal  when  a 
parabolic  pulse  reflects  off  of  a  tilted  surface. 


Ptiltit) 


c  3F  /  2R\ 

\/^isji{(j))Wz^Pw^  \  c  J 

3E  /  _  2R\ 

V^Ap^Wtiit  V  c  J 


(4.50) 


where  Eq.  4.48  was  used.  This  equation  is  rather  complicated.  It  would  be  difficult 
to  derive  CRLBs  or  parameter  estimation  techniques  given  the  signal  in  Eq.  4.50. 
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Parabolic  and  Gaussian  Pulses  Convolution  of  Parabola  and  Gaussian 


Figure  4.7:  (a)  Parabolic  pulses  with  various  widths  and  Gaussian  pulse  with  w  =  2. 

(b)  Convolution  of  parabolic  and  Gaussian  pulses. 

Therefore,  the  remaining  work  in  this  chapter  will  focus  on  the  Gaussian  pulse  model 
from  Sec.  4.2.4  instead  of  the  parabolic  pulse  model. 

4-2.6  Range  Precision  and  the  CRLB.  The  signals  in  Eqs.  4.34,  4.43, 
and  4.50  describe  the  amount  of  optical  power  incident  on  the  LADAR  system’s 
photodetector.  The  photodetector,  typically  a  photodiode  or  photomultiplier  tube, 
converts  the  optical  power  into  an  electrical  signal.  The  relationship  between  the 
optical  power  and  the  photodetector’s  output  current  i(t)  is  given  by  Eq.  2.26 

=  (4.51) 

where  G  is  the  detector’s  gain,  q  is  the  electron  charge,  r]  is  the  fraction  of  photons 
converted  to  electrons  or  quantum  efficiency,  h  is  Planck’s  constant,  and  u  is  the 
laser  optical  frequency.  This  relationship  between  the  received  optical  power  and  the 
photocurrent  is  often  written  in  terms  of  the  detector’s  responsivity  (see  Eq.  2.27). 
Responsivity  has  units  of  charge  per  energy,  or  equivalently,  current  per  power. 
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The  energy  in  the  photocnrrent  signal  (denoted  Eg),  normalized  to  nnit  resis¬ 
tance,  is 

OO  OO 

E,  =  I  e-(t)dt  =  I  PL{t)dt.  (4.52) 

—  OO  — OO 

Notice  that  the  optical  energy  E  is  the  integral  over  Ptmit)  (see  Eq.  4.37)  while  the 
electrical  energy  Ee  is  proportional  to  the  integral  over  Puit{t)  squared.  Thus,  there 
is  not  a  straightforward  general  relationship  between  the  received  optical  energy  E 
and  the  signal  electrical  energy  E^.. 

The  signal  that  is  recorded  by  the  LADAR  system  is  corrupted  by  noise  sources. 
Noise  sources  include  laser  speckle,  background  light  entering  the  LADAR’s  receiver 
aperture,  dark  current  in  the  photodetector,  and  electronic  noise  in  the  analog-to- 
digital  converter  (ADC).  Laser  speckle  and  background  light  are  usually  modeled 
as  Poisson  processes  [21,31].  However,  as  the  mean  of  a  Poisson  random  variable 
increases,  its  distribution  approaches  the  shape  of  a  Gaussian  distribution.  Thus, 
Gaussian  noise  is  sometimes  an  adequate  approximation  for  the  noise  in  LADAR 
systems.  In  this  section,  it  is  assumed  that  the  noise  in  the  received  signal  is  additive 
white  Gaussian. 

The  noisy  received  signal,  denoted  r(t),  is 

r{t)  =  i{t)  +  w{t)  (4.53) 

where  w{t)  is  white  Gaussian  noise  with  a  power  spectral  density  (PSD)  that  has 
height  Nq.  This  continuous-time  signal  is  sampled  at  rate  fs-  The  variance  of  the 
Gaussian  noise  is  denoted  cx^.  If  an  ideal  anti-aliasing  filter  is  used  before  the  analog- 
to-digital  conversion,  then  the  per-sample  noise  variance  is 

a'  =  fsNo.  (4.54) 
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The  sampled  signal  is  denoted  r[n\. 


r[n]  =  r  {n/ fs)  =  i[n]  +  w[n]. 


(4.55) 


Becanse  the  noise  is  additive  white  Ganssian,  the  distribntion  of  the  samples  is 


r[n]  ~  . 


(4.56) 


The  anti-aliasing  hlter  causes  all  of  the  noise  realizations  to  be  independent. 


E  [t(;[m]ia[n]]  =  a'^6[m  —  n] 


(4.57) 


where 

6  [m] 

is  the  Kronecker  delta  function. 


1  if  m  =  0 
0  if  m  7^  0 


(4.58) 


Given  a  signal  mixed  with  AWGN,  the  Gramer-Rao  lower  bound  for  estimates 
of  range  is  [23] 


Var 


R 


> 


c^No 


OO  r  n  2 

dPtiitjt) 


4/ 


dt 


(4.59) 


dt 


where  R  is  any  unbiased  estimate  of  target  range.  Range  estimation  usually  involves 
multiple  steps.  Gombinations  of  matched  filters  and  peak-fitting  algorithms  are  com¬ 
mon.  One  of  the  reasons  that  Gramer-Rao  lower  bounds  are  valuable  is  that  they  are 
valid  regardless  of  the  exact  algorithm  used  to  produce  the  range  estimate.  As  long 
as  the  estimate  produced  by  the  algorithm  is  unbiased,  the  GRLB  is  valid. 


The  signal  model  from  Eq.  4.43  differs  from  the  model  used  in  [23]  in  that  the 
width  of  the  received  signal  is  unknown.  However,  the  Gaussian  pulse  with  AGWN 
model  was  also  studied  in  [16].  In  that  reference,  it  was  proved  addition  of  an  unknown 
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width  parameter  to  the  model  does  not  change  the  range  CRLB.  Thus,  the  CRLB 
from  Eq.  4.59  is  valid  in  this  case. 

Using  Fourier  transform  properties,  the  integral  in  Eq.  4.59  can  be  expressed 

00^00  OO 

dt=  j  \27rfT{Puit{t)}fdf  =  471^  J  f\P{Puitmfdf  (4.60) 

—  OO  — OO  — OO 


r  {dPuitjt) 

J  _  dt 


where  P  {}  denotes  the  Fourier  transform  operation.  The  mean  squared  bandwidth 
of  the  signal,  which  is  denoted  a'j,  is 


OO 

/  p\F{PuHm"<i! 

2  _  -OO 
^  f  OO 

/  \P{P,m{t)}\‘^df 

—  OO 


(4.61) 


Using  this  result,  the  CRLB  can  be  written 


Var 


R 


> 


SNR  \4Traf 


(4.62) 


where  the  SNR  is 

E 

SNR=^.  (4.63) 

No 

It  is  the  ratio  of  the  signal’s  electrical  energy  to  the  height  of  the  AWGN  PSD 

(Afo). 


4-2.7  Gaussian  Pulse  CRLB.  Given  the  Gaussian  pulse  from  Eq.  4.38,  the 
noiseless  electrical  signal  is 


tit)  =  , _ —  exp 

^Wt 


2vji 


2R\ 


The  electrical  energy  in  this  signal  is  (see  Eq.  4.52) 


Ee=  ^  . 

2^wt 


(4.64) 


(4.65) 
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Eq.  4.43  is  used  to  calculate  the  mean  squared  bandwidth.  The  Gaussian  pulse  mean 


square  bandwidth  is 

2  _  1 


(4.66) 


Therefore,  the  CRLB  for  target  range  measurements  made  with  a  Gaussian  pulse  is 


Var 


> 


1  c^w'f 


(4.67) 


Examples  of  Gaussian  pulse  range  GRLBs  are  shown  in  Fig.  4.8.  In  this  hgure, 
a  2.5  cm  beam  intercepts  a  target  10  km  away  from  the  LADAR.  A  1.6  /rm  wavelength 
is  assumed.  In  Fig.  4.8a,  the  transmitted  pulse  FWHM  is  1  ns.  The  increase  in  the 
GRLB  for  a  15  degree  tilt  is  small.  However,  once  the  tilt  is  45  degree,  the  GRLB, 
which  in  standard  deviation  is  proportional  to  the  received  pulse  width  Wt,  has  nearly 
doubled.  For  a  75  degree  tilt,  the  GRLB  is  over  six  times  greater  than  if  the  target 
was  perpendicular  to  the  LOS.  In  Fig.  4.8b,  the  transmitted  pulse  FWHM  is  2  ns. 
This  wider  pulse  is  affected  less  by  target  tilt  than  the  1  ns  pulse.  For  example,  the 
15  and  30  degree  tilts  are  almost  negligible  in  this  case.  However,  for  equal  SNRs, 
the  wider  2  ns  pulse  had  a  worse  range  precision  in  the  absence  of  target  surface  tilt. 
Thus,  given  a  fixed  SNR,  it  is  best  to  use  the  shortest  pulse  possible  even  though 
shorter  pulses  are  more  suspectable  to  distortion  because  of  target  tilt. 

It  was  shown  in  Fig.  4.1  that  if  the  noise  in  a  signal  is  AWGN  that  the  GRLB  can 
be  achieved  unless  the  background  noise  dominates.  However,  the  true  received  pulse 
width  must  be  known  in  order  to  reach  the  bound.  If  the  pulse  width  is  unknown, 
then  the  filter  used  for  range  estimation  will  not  be  designed  properly  and  the  range 
estimates  will  less  precise  than  the  GRLB  predicts.  This  issue  was  studied  for  Poisson 
noise  and  is  illustrated  in  Fig.  3.13. 

The  detrimental  effect  of  target  tilt  on  range  precision  can  be  mitigated  by 
increasing  the  transmitted  beam  size.  Larger  transmitted  beams  make  it  possible  to 
focus  the  beam  into  a  smaller  area  at  range.  Thus,  there  is  less  spreading  of  the 
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Range  CRLBs  for  Tilted  Surfaces 


Range  CRLBs  for  Tilted  Surfaces 


(a) 


(b) 


Figure  4.8:  (a)  Range  CRLBs  for  1  ns  FWHM  pulse, 

(b)  Range  CRLBs  for  target  tilted  45  degrees. 


received  pulse.  The  effect  is  illustrated  in  Fig.  4.9a.  The  target  is  assumed  to  have  a 
45  degree  tilt.  Doubling  the  beam  size  from  2.5  cm  to  5  cm  decreases  the  CRLB  by 
about  one-third.  However,  greater  increases  in  beam  size  yield  a  diminishing  return 
since  the  transmitted  pulse’s  1  ns  FWHM  is  not  decreased.  On  the  other  hand, 
received  optical  energy  is  proportional  to  receiver  aperture  area.  If  the  increase  in 
beam  size  is  accompanied  by  an  increase  in  receiver  apertnre  diameter,  then  the  SNR 
will  increase.  It  is  necessary  to  focus  the  beam  to  gain  this  benefit.  Increasing  the 
size  of  a  collimated  beam  will  actnally  decrease  the  range  precision  since  propagated 
beam  size  for  collimated  beams  increases  monotonically  with  increases  in  transmitted 
beam  size. 

In  Fig.  4.9b,  the  effect  of  varying  the  transmitted  pulse’s  duration  is  stndied.  In 
the  absence  of  target  tilt,  the  range  standard  deviation  is  proportional  to  the  FHWM. 
However,  target  tilt  makes  that  relationship  nonlinear.  Given  a  target  tilted  at  45 
degrees  and  a  2.5  cm  beam  propagating  10  km,  there  is  a  good  increase  in  precision 
when  the  FWHM  drops  from  10  to  3  ns  and  again  when  it  drops  from  3  to  1  ns. 
However,  decreasing  the  FHWM  from  1  ns  to  0.3  ns  does  little  to  improve  the  range 
precision  for  the  45  degree  target.  Withont  tilt,  the  range  standard  deviation  wonld 
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Range  CRLBs  for  Various  Beam  Sizes  Range  CRLBs  for  Various  Pulse  Widths 


(a)  (b) 


Figure  4.9:  (a)  Range  CRLBs  for  1  ns  FWHM  pulse. 

(b)  Range  CRLBs  for  2  ns  FWHM  pulse. 

have  dropped  by  a  factor  of  0.3.  The  performance  difference  between  0.3  and  0.1  ns 
FHWM  pulses  is  almost  not  noticeable  since  in  that  case  the  width  of  the  received 
pulse  is  dominated  by  the  effect  of  target  depth. 

4  •  3  Chapter  Summary 

The  range  precision  of  LADAR  systems  was  discussed  in  this  chapter.  Different 
signal  and  noise  models  were  used  to  analyze  range  precision.  In  Sec.  4.1,  the  range 
CRLB  derived  in  Chapter  III  for  a  truncated  parabolic  pulse  in  Poisson  noise  was 
compared  to  the  bound  for  the  same  pulse  shape  and  Gaussian  noise.  It  was  shown 
that  the  Poisson-noise  CRLB  was  higher  (worse  range  precision)  than  the  Gaussian 
noise  bound.  The  effect  of  pulse  averaging  was  also  discussed  in  Sec.  4.1.  It  was  shown 
that  if  a  hxed  amount  of  laser  energy  is  divided  equally  into  multiple  laser  pulses,  then 
the  range  precision  for  the  accumulated  or  averaged  signal  drops  as  the  number  of 
pulses  increases.  Thus,  transmitting  laser  energy  in  as  few  pulses  as  possible  will 
maximize  a  LADAR’s  range  precision. 

In  Sec.  4.2,  it  was  noted  that  the  reflection  process  can  elongate  a  laser  pulse.  It 
was  also  shown  that  pulse  elongation  decreases  range  precision.  Given  a  pulse  with  a 
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Gaussian  shape  mixed  with  AWGN  noise,  the  range  GRLB  was  calculated  in  terms  of 
target  surface  orientation  and  other  LADAR  system  parameters  such  as  transmitter 
aperture  size,  range  to  target,  and  laser  beam  focus  distance.  Given  an  aperture  that 
is  a  few  inches  in  diameter  and  a  laser  pulse  that  is  at  least  a  nanosecond  in  duration, 
the  effect  of  a  target  surface  tilt  on  range  precision  is  negligible  at  ranges  less  than  a 
few  kilometers.  However,  airborne  LADAR  systems  can  encounter  scenarios  where  a 
target  surface  tilt  can  produce  a  signihcant  drop  in  range  precision.  The  amount  of 
increase  in  the  range  GRLB  due  to  target  tilt  was  quantihed  in  Sec.  4.2. 
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V.  Gain  Variation  in  APD  Arrays 

Some  modern  LADAR  systems  make  use  of  compact  arrays  of  avalanche  photodiodes 
to  produce  high-resolution  three-dimensional  images  of  scenes.  The  arrays  in  those 
systems  are  often  very  compact.  Systems  with  over  ten  thousand  detector  elements 
have  been  demonstrated.  The  sheer  number  of  pixels  present  in  modern  photodetector 
arrays  can  eliminate  the  need  for  scanning  equipment  in  LADAR  systems.  However, 
the  miniaturization  process  used  to  make  the  compact  detector  arrays  can  introduce 
undesirable  effects.  An  unwanted  property  that  has  been  observed  in  arrays  of  APDs 
is  gain  variation. 

Gain  variation  refers  to  fluctuation  in  the  relationship  between  the  number  of 
electrons  detected  by  the  photodiode  and  the  output  current  of  that  device.  Ideally, 
this  relationship  would  be  fixed.  In  practice,  it  has  been  shown  that  gain  can  fluctuate 
over  time.  These  fluctuations  effect  the  entire  array.  The  fluctuation  can  distort  the 
shape  of  the  received  signal  and  can  even  create  false  signals  on  detector  pixels  that 
did  not  actually  observe  a  laser  signal. 

In  this  chapter,  mathematical  models  are  proposed  to  represent  gain  variation  in 
arrays  of  photodiodes.  These  models  of  signals  observed  by  LADAR  detector  arrays 
are  used  to  derive  numeric  methods  of  gain  estimation.  The  numeric  methods  are 
tested  using  Monte  Carlo  simulations  and  LADAR  data. 

In  Sec.  5.1,  three-dimensional  LADAR  data  collected  during  an  experiment 
using  a  large  target  panel  is  displayed  and  the  effect  of  gain  variation  on  range- 
resolved  signals  is  illustrated.  In  Sec.  5.2,  a  numeric  method  of  gain  estimation  based 
on  the  Gauss-Newton  iteration  is  derived  and  tested  on  simulated  and  real  LADAR 
data.  Finally,  in  Sec.  5.3,  an  alternative  method  of  gain  estimation  using  target-free 
LADAR  data  is  derived  and  tested  on  LADAR  data. 

5.1  FLASH  LADAR  Data 

Gain  variation  has  been  observed  in  the  Advanced  Scientific  Concepts  (ASC) 
LADAR  system  operated  by  AFRL  Sensors  Directorate.  The  likely  cause  of  the  gain 
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variation  in  this  LADAR  is  the  fact  that  the  APD  array  was  built  around  one  common 
voltage  regulator  [37].  The  voltage  regulator  attempts  to  keep  the  bias  voltage  Vbias 
constant.  However,  the  bias  voltage  appears  to  fluctuate  over  time.  Since  the  gain  in 
a  photodiode  is  a  function  of  the  bias  voltage  (Eq.  2.50),  the  gain  varies  with  changes 
in  that  voltage. 

The  ASC  LADAR  is  a  flash  system.  Flash  LADAR  systems  use  a  single  laser 
pulse  to  image  an  entire  scene.  The  imaging  is  accomplished  using  an  array  with 
many  pixels.  The  ASC  LADAR  uses  an  128  x  128  array  of  InGaAs  photodiodes.  The 
laser  wavelength  of  this  system  is  1.54/im.  The  detector  array  is  bonded  to  a  readout 
integrated  circuit  (ROIC).  This  ROIC  samples  the  received  LADAR  signals  at  a  rate 
of  420  MHz,  corresponding  to  a  round-trip  range  sampling  resolution  of  0.36  m. 

The  laser  transmitter  produces  1.7  mJ  pulses  with  a  full-width  half-maximum 
of  4.7  ns.  These  pulses  are  diffused  over  an  angle  of  1.5  degrees.  The  LADAR  data 
was  collected  using  a  500  mm  focal  length  lens  with  a  80  mm  aperture.  The  laser  is 
eye-safe  at  ranges  more  than  a  few  inches  from  the  aperture. 

The  APDs  in  the  detector  array  are  separated  by  100  fim.  They  are  approxi¬ 
mately  10  fim  in  width.  An  array  of  lenslets  is  used  to  improve  the  effective  hll  factor 
of  the  detector  array.  With  the  500  mm  lens,  the  instantaneous  field  of  view  (IFOV) 
of  the  pixels  is  0.2  milliradians. 

The  ASC  LADAR  stores  twenty  samples  per  pixel.  Therefore,  the  sampling 
depth  for  the  system  is  only  7.2  m.  The  entire  system  operates  at  a  20  Hz  rate,  which 
provides  three-dimensional  image  information  at  a  rate  compatible  with  typical  video 
displays.  Descriptions  of  LADAR  systems  built  by  ASC  and  examples  of  their  imaging 
capabilities  are  shown  in  [47,48]. 

A  three-dimensional  image  of  a  target  panel  is  shown  in  Fig.  5.1.  The  data  was 
collected  at  Wright-Patterson  Air  Force  Base.  This  image  was  collected  using  the  Air 
Force  Research  Laboratory’s  ASC  LADAR  system.  The  LADAR  system  was  located 
in  a  tower  and  was  23  m  above  the  ground.  The  range  to  the  target  panel  was  70 
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Figure  5.1: 


m.  Thus,  the  downlook  angle  was  20  degrees.  The  target  panel  was  made  of  plywood 
and  was  not  painted  or  coated.  That  panel  was  155  cm  high  and  had  a  similar  width. 

In  Fig.  5.1,  the  intensity  from  the  tenth  of  the  twenty  range  slices  is  shown.  The 
target  panel  was  supported  manually.  During  the  data  collections,  the  target  panel 
orientation  with  respect  to  the  LOS  was  varied.  There  is  signihcant  fluctuation  in  the 
individual  signals,  so  they  are  averaged  over  multiple  image  frames.  In  this  case,  108 
signals  were  averaged. 

Range-resolved  LADAR  data  collected  by  the  ASC  system  is  illustrated  in  Fig. 
5.2.  In  Fig.  5.2a,  four  signals  from  the  target  panel  are  illustrated.  In  this  experiment, 
the  target  panel  was  tilted.  Data  with  lower  row  numbers  corresponded  to  longer 
ranges  to  the  panel.  For  all  four  signals,  the  hrst  hve  or  six  samples  are  noise.  In  row 
one-hundred,  the  target  appears  in  the  six  or  seventh  sample.  The  row  forty  target 
was  not  reached  until  sample  ten.  Near  samples  eight  and  nine,  a  drop  in  gain  is 
visible  in  the  row  forty  and  row  sixty  data  before  the  target  panel  is  observed  in  those 


126 


Range-Resolved  LADAR  Signals 


Range-Resolved  LADAR  Signals 


Figure  5.2:  (a)  LADAR  signals  from  target  panel, 

(b)  Background  noise  in  LADAR  image. 


pixels.  The  variation  in  gain  also  altered  the  observed  width  of  the  laser  pulse.  The 
pulse  is  much  wider  in  rows  eighty  and  one-hundred  than  in  rows  twenty  and  forty. 

Fig.  5.2b  illustrates  data  from  pixels  where  no  target  is  present.  In  these 
pixels,  the  first  seven  samples  contain  only  noise  with  mean  level  around  five  hundred. 
However,  in  samples  eight  and  higher,  the  gain  has  dropped  and  the  mean  is  around 
two  hundred.  This  data  is  typical  of  what  was  observed  throughout  the  experiment. 
A  laser  signal  or  signals  affected  the  gain  on  all  pixels  in  the  array.  The  gain  variation 
was  observed  whether  or  not  there  was  actually  a  laser  signal  in  the  pixel. 


5.2  Gauss-Newton  Method  for  Gain  Estimation 

In  this  section,  the  Gauss-Newton  method  of  parameter  estimation  is  used  to 
solve  for  the  unknown  gain  in  a  LADAR  detector  array.  This  method  requires  that 
the  probability  density  function  of  the  data  is  known.  Therefore,  a  signal  model 
is  developed  to  approximate  the  data  observed  by  a  LADAR  system.  This  model 
assumes  white  Gaussian  noise,  but  does  not  assume  knowledge  of  the  signal  amplitudes 
or  the  target  ranges.  Those  parameters  are  estimated  along  with  the  unknown  gain 
values  by  the  Gauss-Newton  iteration. 
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5.2.1  Detector  Array  Signal  Model.  The  model  used  to  represent  signals 
observed  by  the  detector  array  is  based  on  the  LADAR  signal  model  from  [4],  How¬ 
ever,  some  modifications  have  been  made  to  that  equation.  The  temporal  pulse  shape 
is  assumed  to  be  a  Gaussian  function,  rather  than  the  truncated  parabola  from  Eq. 
2.29.  Also,  it  is  assumed  that  the  noise  in  all  samples  is  modeled  by  independent  and 
identically  distributed  (HD)  Gaussian  random  variables.  Some  detector  signal  models 
used  earlier  in  this  dissertation  assumed  Poisson  noise. 

The  detector  array  signal  model  is 


Ik,m 

=  GklAra  exp 


{k  —  R„ 


+  -Bn 


(5.1) 


where 


ik,m  =  Ara  exp 


{k  —  R„ 


-Bn 


(5.2) 


is  the  signal  before  the  variable  gain  is  applied.  In  Eq.  5.1,  k  denotes  the  temporal 
sample  index  and  m  denotes  the  detector  pixel  (spatial)  index.  Note  that  in  this 
equation,  the  fc,  i?™.,  and  w  have  units  of  samples.  Each  detector  records  K  samples 
and  there  are  M  detectors  in  the  array.  The  variable  Am  is  the  peak  amplitude  of  the 
received  signal  before  gain.  The  range  to  the  target,  in  samples,  is  denoted  Rm-  The 
bias  level  on  the  m-th  detector  element  is  Em¬ 


it  is  assumed  that  the  gain  values  (the  G^’s)  are  unknown.  It  is  also  assumed 
that  the  signal  amplitude  and  range  values  (the  Am’s  and  the  -Rm’s)  are  unknowns. 
However,  the  biases  on  each  detector’s  signal  (the  -B^’s)  are  assumed  to  be  known.  In 
a  real  LADAR  system,  the  bias  is  not  known  and  must  be  estimated  from  the  data. 
However,  noise  data,  which  can  be  used  to  measure  the  bias,  is  normally  abundant. 
The  large  amount  of  available  noise  data  can  be  used  to  produce  precise  bias  estimates 
for  each  detector  in  the  array. 
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The  pulse  width  w  is  assumed  known.  There  are  LADAR  imaging  scenarios 
where  this  assumption  is  not  realistic.  However,  during  the  data  collection  with  the 
ASC  system,  the  range  to  the  target  was  short.  Also,  the  system’s  instantaneous  field 
of  view  is  narrow.  Under  these  circumstances  (short  range  and  narrow  IFOV),  target 
depth  in  individual  pixels  is  normally  shallow  and  the  pulse  width  does  not  stretch 
significantly  because  of  reflection. 

The  k-th  sample  recorded  by  the  m-th  detector  in  the  array  is  the  random 
variable  Dk,m-  This  random  variable  has  mean  Ik^m  and  additive  noise.  It  has  been 
assumed  that  the  noise  is  Gaussian.  The  variance  of  each  sample  is  a^. 


Dk,m  ~  A/"  ^  ]  • 


(5.3) 


The  probability  density  function  is 


PDk,m  {dk,m) 


72 


:  exp 


TTCr^ 


2ct2 


(5.4) 


where  dk^m  is  a  realization  of  random  variable  Since  the  noise  is  IID,  the  joint 

distribution  of  the  data  is  the  product  of  the  PDFs  of  the  dfc,m’s.  Let  D  denote  the  set 
of  random  variables  and  d  denote  the  realization  of  the  set  of  those  random  variables. 
The  joint  PDF  of  D  is 


K  M 


po  (d)  =  n  n  (*',m') 
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Let  g  denote  the  set  of  K  unknown  gain  values,  let  a  be  the  set  of  M  unknown 
amplitude  values,  and  let  r  be  the  set  of  M  range  values.  The  log-likelihood  function 
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for  random  data  set  d  is 


I  (g,  a,  r|d)  =  log  {2TTa^)  -  ^  ^  ^  f  ■  (5.6) 

k'=l  m'=l 

This  log-likelihood  function  is  relatively  complicated.  Its  complexity  will  likely  prevent 
efficient  calculation  of  maximum  likelihood  estimates.  However,  in  the  next  section, 
it  is  shown  that  it  can  be  used  in  combination  with  the  Gauss-Newton  method  to 
derive  iterative  estimates  of  the  unknown  parameters. 


5.2.2  Iterative  Gain  Estimation.  The  unknown  parameters  in  the  detector 
array  signal  model  from  Eq.  5.1  can  be  estimated  using  the  Gauss-Newton  method. 
The  Gauss-Newton  method  is  an  iterative  technique  for  parameter  estimation  that  is 
based  on  the  log-likelihood  function.  The  Gauss-Newton  method  is  described  in  Fun¬ 
damentals  of  Statistical  Signal  Processing:  Estimation  Theory  by  Kay  [23] .  Given  a 
log-likelihood  function  with  unknown  parameter  vector  0,  the  Gauss-Newton  method 
relies  on  an  initial  estimate  of  the  unknown  parameters,  Oq.  This  estimate  is  not 
expected  to  be  an  optimal  measurement  of  the  unknown  parameters.  However,  the 
Gauss-Newton  method  requires  an  precise  initial  estimate.  If  the  initial  estimate  is  not 
precise,  the  iteration  could  converge  to  a  local  maximum  as  opposed  to  the  desired 
global  maximum  of  the  likelihood  function.  Given  estimate  0„,  the  Gauss-Newton 
iterative  update  is 


^n+l  — 


-1 


A(0n)  A(0„). 


(6.7) 


The  matrix  A  and  the  vector  A  in  Eq.  5.7  are  populated  by  derivatives  of  the  log- 
likelihood  function.  The  vector  A  contains  the  hrst  derivatives  of  the  log-likelihood 
function.  The  Tth  element  of  A  is  the  derivative  of  I  (6)  with  respect  to  the  Gth 
unknown  parameter. 


[A]. 


dl  (0) 


(5.8) 


where  6i  is  the  i-th  element  of  0.  The  elements  of  matrix  A  are  the  second  derivatives 
of  the  log-likelihood  function.  The  element  in  the  Tth  row  and  j-th  column  of  matrix 
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(5.9) 


A  is 


dH  (6) 

dOidOj 


The  likelihood  function  is  parameterized  by  6.  That  parameter  vector  is  being 
updated  iteratively.  Within  the  iteration,  A  and  A  are  evaluated  using  the  unknown 
parameter  vector  estimate  from  the  previous  step.  This  explains  the  notation  A(0„) 
and  A(0n)  in  Eq-  5.7. 

The  hrst  derivatives  of  the  log-likelihood  function  are 
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(5.10) 


The  second  derivative  with  respect  to  6i  is 
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(5.11) 


These  second  derivatives  are  the  diagonal  elements  of  A.  Matrix  A  also  contains 
second  derivatives  with  respect  to  9i  and  9j.  These  derivatives  are  the  off-diagonal 
elements  of  A.  These  derivatives  are 
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(5.12) 


These  log-likelihood  derivative  expressions  with  respect  to  general  unknown  param¬ 
eters  9i  and  9j,  can  be  used  to  find  elements  of  A  and  A  specific  to  the  likelihood 
function  associated  with  the  joint  density  from  Eq.  5.5. 

The  next  step  in  the  derivation  of  the  Gauss-Newton  iteration  is  evaluation  of 
the  derivatives  of  the  signal  from  Eq.  5.1.  The  derivative  of  the  signal  with  respect 
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to  the  fc-the  gain  value  is 


(5.13) 


dh 


k,m 


dG, 


The  signal  Ik^m  does  not  depend  on  the  gain  on  samples  other  than  the  fc-th  sample. 
Therefore, 


dL 


k,m 


dGi 


=  0  when  k  ^  1. 


(5.14) 


The  signal  derivatives  with  respect  to  the  amplitude  and  range  for  the  m-th  detector 
are 
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and 
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The  amplitude  and  range  of  signal  h^^^  do  not  depend  on  range  and  amplitude  values 
of  signals  on  other  detectors.  Therefore, 


dh 
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=  0  when  m  ^  n. 


(5.17) 


Second  derivatives  of  are  also  needed  to  evaluate  the  terms  in  the  Gauss- 
Newton  iteration.  Second  derivatives  with  respect  to  gain  and  amplitude  are  zero. 
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The  range  second  derivative  is  not  zero. 
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Terms  dependent  on  derivatives  of  different  gain  valnes  are  zero 


d^L 


k^m 


OGkdGi 


=  0  when  k  ^  1. 


(5.20) 


Similarly,  derivatives  dependent  of  different  amplitude  or  different  range  values  are 
zero. 
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=  0  when  m  ^  n. 
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The  derivatives  corresponding  to  gain  and  amplitude  values  are 
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(5.22) 


The  gain  and  range  derivatives  are 
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The  final  derivative  needed  is 


-  —  (A;  -  i?  )  exp 

dA^m^ 


(fc  -  RmY 

2w‘^ 
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The  first  derivatives  of  /(g,  a,  r|d)  are  needed  to  evaluate  the  elements  of  A. 
Eqs.  5.8  and  5.10  are  combined  to  find  the  derivative  with  respect  to  unknown  gain 
terms  Gk- 
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where  Eqs.  5.13  and  5.14  were  used.  The  other  elements  of  A  are 
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and 
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(5.27) 


The  hrst  and  second  derivatives  of  I  (g,  a,  r|d)  are  needed  to  evaluate  the  ele¬ 
ments  of  A.  The  diagonal  elements  of  A  are  evaluated  hrst,  beginning  with  the  gain 
terms.  Using  Eq.  5.11 
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The  terms  associated  with  the  unknown  signal  amplitudes  (the  A^’s)  are 


(5.28) 


d^l  (g,a,r|d) 

dAl 


K 


CT" 


fc/exp 


fc'=i 


jk'  -  RmY 

2w‘^ 


(5.29) 


The  range  terms  are 
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Now  the  off-diagonal  elements  of  A  are  calcnlated.  Terms  corresponding  to 


derivatives  to  different  gain  valnes  are  zero. 


d‘^1  (g,a,r|d) 

dGkdGi 


=  0  when  k  ^  1. 


(5.31) 


Several  other  off-diagonal  elements  of  A  are  zero. 
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However,  terms  derived  using  derivatives  of  gain  and  amplitude  are  not  zero.  These 
terms  are 
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The  gain  and  range  terms  are 


(5.33) 
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The  hnal  terms  from  A  are 
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All  terms  necessary  to  evaluate  A  and  A  have  been  derived.  The  elements  of 
A  are  shown  in  Eqs.  5.25-5.27.  The  elements  on  the  diagonal  of  A  are  found  in 
Eqs.  5.28-5.30.  The  off-diagonal  elements  A  are  listed  in  Eqs.  5.31-5.35.  In  the  next 
section,  the  Gauss-Newton  iterative  solution  to  the  array  signal  model  is  tested  via 
simulation. 
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5.2.3  Simulations.  In  this  section,  the  Gauss-Newton  iterative  gain  esti¬ 
mation  process  derived  in  Sec.  5.2.2  is  tested  using  Monte  Carlo  simulations.  The 
simulated  signals  are  illustrated  in  Fig.  5.3.  Four  signals  were  generated  (M  =  4). 
Each  of  those  signals  had  thirty-two  samples  {K  =  32).  The  range  and  amplitudes  of 
the  signals  were  generated  randomly.  The  range  values,  measured  in  samples,  were 
between  eight  and  twenty-four.  The  amplitude  values  were  between  ten  and  twenty. 

The  width  of  the  Gaussian  pulse  was  set  at  three  samples  {w  =  3).  This  value 
was  assumed  to  be  known  in  the  Gauss-Newton  iteration.  The  bias  level  on  each 
signal  was  generated  randomly.  The  bias  values  were  between  five  and  ten.  These 
values  were  assumed  to  be  known  for  the  iteration.  In  a  real  LADAR  system,  the 
background  or  bias  level  must  be  estimated.  However,  noise  data  is  abundant  and 
precise  background  level  measurements  are  usually  available. 

The  additive  IID  Gaussian  noise  in  each  sample  had  unit  variance  (cr^  =  1). 
Knowledge  of  this  value  is  not  needed  by  the  Gauss-Newton  iteration.  All  elements  of 
vector  A  and  A  are  inversely  proportional  to  cx^.  Therefore  A  and  A  are  also  inversely 
proportional  to  a^.  The  Gauss-Newton  iteration  (Eq.  5.7)  contains  the  term  A~^A. 
The  variance  in  that  term  cancels  since  A“^  is  proportional  to  Thus,  the  variance 
is  not  needed  in  the  Gauss-Newton  iteration  equations. 

A  series  of  gain  values,  common  to  all  four  signals,  was  generated.  The  sequence 
of  gain  values  in  the  Monte  Carlo  was  a  sinusoid  with  values  between  0.8  and  1.2. 
The  phase  of  this  sinusoid  was  generated  randomly  and  the  frequency  was  selected  so 
that  there  were  about  two  complete  oscillations  over  the  thirty-two  samples. 

In  Fig.  5.3a,  the  simulated  signals  before  the  gain  is  applied  are  illustrated  by 
solid  lines.  The  signals  after  the  multiplicative  gain  values  are  applied  are  shown  as 
broken  lines.  The  signals  before  and  after  gain  are  related  according  to  Eq.  5.1.  There 
is  also  a  random  realization  of  simulated  LADAR  data  in  Fig.  5.3b.  The  signals  with 
gain  (broken  lines)  are  shown  for  comparison.  The  unit  variance  additive  noise  is 
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Detector  Array  Signal  Model 


(a) 


Detector  Array  Signal  Model 


Sample 

(b) 


Figure  5.3:  (a)  Simulated  signals  before  and  after  variable  gain, 

(b)  Noisy  data  with  variance  =  1. 


distributed  according  to  Eq.  5.5.  The  noise  is  noticeable,  but  is  not  very  powerful 
compared  to  the  signal  energy. 

The  simulated  LADAR  data  shown  in  Fig.  5.3  was  processed  using  the  Gauss- 
Newton  iteration  from  Eq.  5.7.  That  equation  required  an  initial  estimate  of  the 
unknown  parameter  vector  9q.  The  actual  order  in  which  the  gain,  amplitude,  and 
range  estimates  are  arranged  in  0  is  arbitrary.  The  convention  adopted  in  this  disser¬ 
tation  is  to  order  the  gain  estimates  hrst,  then  the  amplitude  estimates,  and  hnally 
the  range  estimates. 

The  initial  K  gain  estimates  are  all  assumed  to  be  one. 


1  V  k. 


(5.36) 


The  M  initial  amplitude  estimates  are  determined  from  the  data.  These  estimates 
are  the  peak  value  on  each  detector  minus  the  known  bias  on  that  detector. 
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The  initial  range  estimate,  in  samples,  is  the  argument  k  corresponding  to  the  maxi¬ 
mum  recorded  value  on  the  m-th  detector. 


ra+2K 


arg  max  {d^^rn  '■  k 


1,2,. 
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The  initial  estimates  from  Eqs.  5.36,  5.37,  and  5.38  are  used  to  initialize  the 
Gauss-Newton  iteration.  Those  initial  estimates  are  also  needed  to  evaluate  A  and  A 
in  the  Gauss-Newton  iteration.  The  initial  iteration  is 
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A(0o) 


A(0o). 


(5.39) 


After  this  initial  calculation  is  made,  the  process  is  repeated  using  Eq.  5.7.  In 
theory,  this  Gauss-Newton  iteration  always  converges.  However,  it  is  possible  for 
it  to  converge  to  a  local  maximum  rather  than  the  desired  global  maximum  of  the 
likelihood  function.  In  practice,  numeric  errors  in  the  calculations  can  make  the 
iteration  unstable.  In  that  case,  the  parameter  estimates  may  fail  to  converge  to 
finite  values. 

It  is  straightforward  to  tell  when  the  iteration  is  not  converging.  In  that  case, 
the  magnitudes  of  at  least  some  of  the  estimates  approach  inhnity.  Also,  the  condition 
number  of  matrix  A  becomes  very  large.  The  condition  number  is  defined  to  be  the 
ratio  of  the  largest  to  the  smallest  singular  value  of  a  matrix  [18].  If  the  Gauss-Newton 
iteration  does  not  converge  it  can  be  tried  again  using  different  initial  estimates.  In  a 
LADAR  system  with  multiple  detectors,  it  is  also  possible  to  omit  data  from  certain 
detectors  when  implementing  the  iteration.  It  is  possible  that  omitting  data  from 
some  detectors  could  cause  the  iteration  to  converge  even  though  it  did  not  when 
data  from  all  detectors  was  used. 

Using  the  simulated  LADAR  data  shown  in  Fig.  5.3,  the  Gauss-Newton  derived 
in  Sec.  5.2.2  iteration  was  tested.  The  initial  estimate  of  0  was  calculated  using 
Eqs.  5.36,  5.37,  and  5.38.  Twenty  iterations  of  the  algorithm  were  computed.  The 
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amplitude  and  range  estimates  for  each  of  the  four  simulated  detectors  are  shown  in 
Fig.  5.4. 

Range  estimation  was  straightforward  since  the  initial  estimates  were  precise. 
All  initial  estimates  were  within  two  samples  of  the  true  value  (see  Fig.  5.4a).  The 
algorithm  converged  to  stable  and  precise  range  estimates  with  a  few  steps.  Amplitude 
estimation  was  more  challenging  because  of  the  unknown  gain  values  in  the  simulated 
LADAR  data.  Initial  amplitude  estimates  were  off  by  over  one-quarter  of  the  true 
value.  However,  after  about  ten  steps,  the  Gauss-Newton  iteration  converged  to  stable 
and  relatively  precise  estimates  of  the  true  amplitude  values  (see  Fig.  5.4b). 

The  iteration  also  yielded  estimates  of  the  unknown  gain  values  in  the  simulated 
LADAR  data.  The  hnal  gain  estimates  after  twenty  iterations  of  the  Gauss-Newton 
algorithm  are  shown  in  Fig.  5.5a.  Most  of  the  estimates  are  within  ten  percent  of 
the  true  gain  value  even  though  only  four  detectors  were  simulated.  If  there  had  been 
data  from  more  than  four  detectors,  then  these  gain  estimates  would  have  been  even 
more  precise. 

Using  the  hnal  estimates  of  gain,  amplitude,  and  range,  the  original  LADAR 
signals  can  be  estimated.  These  estimates  are  plotted  in  Fig.  5.5b.  The  true  signals 
are  shown  for  comparison.  There  is  good  agreement  between  the  true  signals  and 
the  estimated  signals  that  are  dehned  by  the  output  parameter  estimates  from  the 
Gauss-Newton  iteration. 

5.2.4  Gain  Estimation  with  LADAR  Data.  The  Gauss-Newton  iteration 
derived  in  Sec.  5.2.2  was  tested  on  LADAR  data  collected  by  the  ASG  camera.  Data 
from  the  target  panel  shown  in  Fig.  5.1  was  used.  The  algorithm  was  seeded  using 
the  initialization  equations  from  Sec.  5.2.3. 

The  Gauss-Newton  iteration  was  tested  on  many  different  combinations  of  pixels 
from  the  target  panel.  The  performance  of  the  algorithm  using  real  LADAR  data  was 
bad.  In  the  vast  majority  of  cases,  the  iteration  diverged  immediately.  In  those  cases, 
estimates  of  range,  amplitude,  and  gain  immediately  trended  toward  unreasonable 
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Figure  5.4:  (a)  Iterative  estimation  of  the  unknown  range  values, 

(b)  Iterative  estimation  of  the  unknown  amplitude  values. 


Detector  Array  Signal  Model 


(a) 


Detector  Array  Signal  Model 


(b) 


Figure  5.5: 


(a)  Estimation  of  unknown  gain  values,  (b)  Estimated  signals. 
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LADAR  Data  Gain  Estimation 
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Figure  5.6:  (a)  Gain  estimates  after  ten  iterations. (b)  Signal  estimates. 


values.  When  the  iteration  diverged,  the  condition  number  of  A  became  very  large 
within  a  few  steps  and  that  matrix  became  impossible  to  accurately  invert. 

In  some  tests  of  the  Gauss-Newton  iteration,  the  gain  estimates  started  to  look 
similar  to  the  shape  of  the  range-resolved  LADAR  data.  One  of  these  cases  is  illus¬ 
trated  in  Fig.  5.6a.  The  four  pixels  of  LADAR  data  used  in  this  case  are  shown  in  Fig. 
5.6b.  The  algorithm  was  allowed  to  iterate  ten  times.  The  gain  estimates  correspond¬ 
ing  to  background  (noise)  data  were  approximately  one.  However,  the  gain  estimates 
when  laser  pulses  were  present  in  the  data  were  several  times  greater  than  one  and 
took  on  the  shape  of  the  pulses.  In  this  case,  the  modeled  LADAR  signals  using  the 
algorithm’s  estimates  were  reasonably  accurate.  However,  the  modeled  signals  did 
not  accurately  capture  the  range  variations  in  the  signals.  The  closest  and  farthest 
LADAR  signals  are  about  two  samples  apart.  But  the  range  estimates  provided  by 
the  Gauss-Newton  iteration  are  all  within  a  fraction  of  a  sample. 

Inspection  of  the  range  and  amplitude  values  during  iteration  shows  a  problem 
in  the  process.  The  iteration  consistently  drove  the  signal  amplitude  estimates  toward 
zero  (see  Fig.  5.7b).  In  some  cases,  the  amplitude  estimates  took  negative  values. 
Negative  amplitudes  are  physically  meaningless.  As  illustrated  in  Fig.  5.6,  the  Gauss- 
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Figure  5.7:  (a)  Iterative  estimation  of  the  unknown  range  values, 

(b)  Iterative  estimation  of  the  unknown  amplitude  values. 


Newton  iteration  preferred  to  model  the  laser  pulses  using  the  gain  values  rather  than 
the  Gaussian  pulse  shape  embedded  in  the  signal  model.  As  the  amplitude  values 
dropped,  the  range  estimates  became  unstable.  The  range  estimates  started  to  vary 
quickly  between  steps  of  the  iteration.  Eventually,  the  range  estimates  always  moved 
outside  the  range  of  meaningful  values  and  the  condition  number  of  A  became  very 
large. 

The  Gauss-Newton  iteration  for  gain  estimation  was  tested  on  many  thousands 
of  combinations  of  pixels  of  LADAR  data.  The  gain  estimates  should  have  observed 
a  drop  over  time.  This  drop  in  gain  is  obvious  in  the  background  data  shown  in 
Fig.  5.2.  However,  good  algorithm  performance  was  never  observed.  The  algorithm 
occasionally  modeled  the  LADAR  signals  accurately.  However,  that  only  happened 
when  the  iteration  used  the  gain  values  to  approximate  the  pulse  shape.  The  Gauss- 
Newton  iteration  algorithm  never  successfully  measured  the  gain  variation  in  the 
LADAR  data. 
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5.3  Gain  Equalization  using  Background  Data 

In  this  section,  a  second  method  of  gain  estimation  is  proposed  and  tested.  This 
method  uses  data  from  target-free  (background  light  and  electronic  noise  only)  pixels 
to  observe  the  gain  variation.  The  observed  variation  in  APD  gain  is  used  to  equalize 
LADAR  data  that  has  targets.  It  is  shown  that  this  method  improves  the  precision 
of  range  estimates  made  with  the  LADAR  system. 

A  modihcation  of  the  photocurrent  from  Eq.  4.51  is  used  in  this  section.  The 
modified  version  includes  a  term  for  the  background  or  noise  power  observed.  It  is 
assumed  that  the  background  power  is  constant  in  the  signal.  The  photocurrent  is 
related  to  the  received  optical  power  by 

^{t)  =  ^[Pr{t)+P,\  (5.40) 

where  G  is  the  detector’s  gain,  q  is  the  charge  of  an  electron,  rj  is  the  detector’s 
quantum  efficiency,  h  is  Planck’s  constant,  and  u  is  the  laser’s  optical  frequency.  The 
term  Pb  is  the  power  from  background  light  in  the  scene  and  other  noise  sources  in 
the  system.  This  power  level  is  measured  before  amplification  (gain).  Background 
power  is  assumed  constant  and  is  amplified  in  the  same  way  that  the  received  laser 
signal  power  is  increased. 

Ideally,  the  gain  of  the  detector  would  be  constant.  However,  it  has  been  shown 
in  this  dissertation  and  in  [37]  that  gain  fluctuates  in  APD  arrays.  Gain  variation 
distorts  the  shape  of  received  laser  pulses.  In  imaging  LADAR  systems  using  focal 
planes,  the  gain  variations  can  affect  the  entire  array.  Background  light  is  amplified 
by  APDs  in  the  same  way  that  received  laser  signals  are.  Because  the  entire  array  is 
affected  by  the  gain  variation,  it  is  possible  for  pixels  with  no  target  present  to  appear 
as  though  they  have  observed  one  because  the  signal  level  is  fluctuating. 

A  relatively  simple  gain  equalization  algorithm  is  proposed  and  tested.  This 
algorithm  uses  data  from  target-free  pixels  to  estimate  the  avalanche  photodiode  gain 
over  time.  It  is  shown  that  the  equalization  algorithm  can  decrease  range  error  in 
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LADAR  data.  It  is  possible  that  this  method  could  be  applied  to  design  of  APD 
arrays.  A  single  or  a  few  pixels  could  be  isolated  and  used  for  measurement  of  gain 
variation.  That  data  would  be  available  to  the  ROIC  for  gain  equalization. 

5.3.1  Gain  Variation  in  Flash  LADAR  Data.  The  same  LADAR  data  that 
was  shown  in  Sec.  5.1  and  processed  in  Sec.  5.2.4  is  used  in  this  section.  All  of  the 
three-dimensional  LADAR  images  are  of  a  large  target  panel  that  was  supported 
manually.  In  image  of  the  panel  is  shown  in  Fig.  5.1  and  examples  of  range  resolved 
LADAR  signals  are  shown  in  Fig.  5.8. 

The  data  displayed  in  Fig.  5.8a  is  from  a  three-dimensional  LADAR  image 
where  the  target  panel  was  tilted  -20  degrees  with  respect  to  the  LOS.  The  tilt  angle 
was  dehned  so  that  a  negative  tilt  corresponded  to  the  top  of  the  target  moving  toward 
the  LADAR.  Thus,  a  positive  tilt  implies  that  the  panel’s  top  has  moved  away  from 
the  sensor.  In  the  case  shown  in  Fig.  5.8a,  the  target  panel  was  oriented  perpendicular 
to  the  ground.  Because  the  LADAR  viewed  the  panel  from  a  tower,  the  top  of  the 
panel  is  the  part  that  is  closest  to  the  system.  Range  resolved  signals  are  shown 
in  Fig.  5.8a.  As  the  row  number  of  these  signals  increases,  the  leading  edge  of  the 
recorded  laser  pulses  moves  to  the  right.  The  trailing  edge  is  also  shifted  to  the  right 
with  increasing  row  number.  However,  the  trailing  edge  does  not  shift  as  quickly  as 
the  leading  edge  does. 

The  leading  edges  of  pulses  reflected  off  of  the  top  and  bottom  of  the  target  panel 
are  separated  by  about  2.5  samples.  However,  the  trailing  edges  are  only  separated  by 
about  1.5  samples.  The  received  laser  pulses  should  have  the  same  shape.  However, 
the  gain  on  those  pulses  from  the  avalanche  photodiode  array  varies  over  time.  Thus, 
the  width  of  the  recorded  signals  decreases  as  the  target  range  increases. 

In  Fig.  5.8b,  the  panel  was  tilted  so  that  its  bottom  is  the  part  closest  to 
the  LADAR.  In  this  case,  the  signals  shift  to  the  left  as  the  row  number  increases. 
The  tilt  was  32  degrees.  Because  of  the  target  panel  tilt,  the  leading  edges  of  the 
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Range-Resolved  LADAR  Signals  Range-Resolved  LADAR  Signals 


Figure  5.8:  (a)  Target  panel  tilt  is  -20  degrees,  (b)  Target  panel  tilt  is  32  degrees. 

recorded  pulses  are  separated  by  about  3.5  samples.  However,  the  trailing  edges  are 
only  separated  by  about  2  samples. 

Inspection  of  the  data  displayed  in  Figs.  5.8a  and  5.8b  shows  that  the  LADAR 
signals  recorded  by  the  APD  array  are  distorted  in  a  predictable  manner.  The  signals 
in  pixels  closest  to  the  LADAR  are  have  the  longest  duration.  Signals  recorded  in 
pixels  farther  from  the  LADAR  are  narrower  regardless  of  the  orientation  of  the  panel 
with  respect  to  the  system.  The  cause  of  the  distortion  in  laser  pulse  shape  is  gain 
variation  in  the  array  of  APDs. 

5.3.2  Gain  Variation  Equalization.  It  is  hypothesized  that  the  gain  of 
the  entire  avalanche  photodiode  array  drops  when  laser  signals  are  detected.  The 
relationship  between  background  photocurrent  level  ih{t)  and  gain  is  given  by  Eq. 
5.40  with  signal  power  Pr{t)  =  0. 

ih{t)  =  (5.41) 

Data  from  target-free  pixels  is  shown  in  Fig.  5.9.  There  is  signihcant  fluctuation 
in  the  background  data.  Thus,  it  is  averaged  over  108  frames  and  over  multiple  (thirty) 
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Figure  5.9:  (a)  Background  data  with  panel  tilted  -20  degrees,  (b)  Background 

data  with  panel  tilted  32  degrees. 


pixels.  The  background  signal  level  is  steady  for  the  hrst  six  samples.  After  that,  the 
background  level  drops  by  a  factor  of  more  than  two.  This  drop  in  gain  occurs  at  the 
point  where  the  target  panel  is  located  (see  Fig.  5.8). 

Averaged  background  data  is  used  to  estimate  the  time-varying  gain.  The  gain 
estimate  G,  which  is  written  here  as  a  function  is  time,  is  calculated  using  Eq.  5.41. 


d{t)  oc  isit) 


(5.42) 


where  isit)  is  the  average  of  background  signals.  This  estimate  of  the  gain  is  used  to 
correct  the  variation  in  the  LADAR  signals  shown  in  Fig.  5.8.  The  data  from  those 
plots,  divided  by  the  average  background  (that  is,  i(t)/iB(t))  is  shown  in  Fig.  5.10. 

The  equalized  data  in  Fig.  5.10a  is  from  images  with  the  panel  tilted  at  -20 
degrees  (the  raw  data  is  shown  in  Fig.  5.8a).  Gain  equalization  removed  most  of 
the  width  variation  in  the  pulses.  However,  the  trailing  edge  of  the  row  112  pulse  is 
out  of  position  compared  to  the  pulses  from  rows  54,  73,  and  93.  But  overall,  the 
gain  equalization  process  has  improved  the  quality  of  the  data  in  spite  of  the  row  112 
outlier. 
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Gain  Equalized  LADAR  Signals 


Figure  5.10:  (a)  Target  panel  tilt  is  -20  degrees,  (b)  Target  panel  tilt  is  32  degrees. 


In  Fig.  5.10b,  the  data  from  panel  tilted  32  degrees  from  LOS  is  shown  after 
equalization.  (The  raw  data  is  shown  in  Fig.  5.8b.)  In  this  case,  equalization  worked 
very  well.  The  gain  equalized  laser  pulses  all  have  the  same  width  regardless  of  their 
position  on  the  panel.  Equivalently,  the  pulse  width  no  longer  depends  on  the  range 
to  the  target. 


5.3.3  Range  Measurements.  In  this  section,  raw  and  equalized  LADAR 
signals  are  used  to  make  measurements  of  range  to  the  target  panel.  The  range 
estimates  were  calculated  by  computing  a  threshold  crossing.  Specihcally,  the  point 
where  the  leading  edge  of  the  pulse  is  one-half  the  maximum  height  of  the  pulse  is 
considered  to  be  the  target  range.  The  range  estimate  was  determined  using  a  linear 
ht  to  the  pulse’s  leading  edge.  The  ht  used  the  points  that  were  just  above  and  below 
the  threshold. 

The  range  estimation  algorithm  is  illustrated  in  Fig.  5.11.  The  peak  value  of  the 
signal  is  located  and  a  parabolic  ht  (show  in  red)  is  calculated  to  the  points  around 
the  peak.  For  range  estimation,  the  threshold  is  set  at  one-half  of  the  peak  ht  value 
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LADAR  Signals  (Column  #65) 


Figure  5.11;  Illustration  of  range  estimation  algorithm. 


(shown  as  a  broken  green  line).  The  range  estimate  is  the  location  where  the  data 
hrst  crosses  the  threshold  (shown  as  a  solid  green  line). 

It  is  common  to  use  a  matched  filter  when  processing  range-resolved  data.  How¬ 
ever,  the  use  of  matched  hlters  assumes  that  the  shape  of  the  received  pulse  is  known. 
Matched  hlters  also  assume  that  the  background  level  is  not  changing.  The  gain  vari¬ 
ation  is  obviously  distorting  the  pulse  shape,  making  the  use  of  a  true  matched  hlter 
impossible. 

It  has  also  been  shown  that  three-dimensional  imaging  LADARs  with  short 
range  sampling  intervals  such  as  this  one  sometimes  fail  to  sample  the  received  signal 
completely  [4].  The  pulses  shown  in  Fig.  5.8  take  up  most  of  the  range  record.  If  the 
trigger  had  been  off  by  a  few  nanoseconds,  these  pulses  would  not  have  been  sampled 
completely.  When  the  pulse  is  sampled  incompletely,  correlation-based  operations 
such  as  matched  hltering  do  not  function  as  intended. 
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Because  of  the  gain  variation,  the  distortion  of  the  received  laser  pulses,  and 
the  fact  that  pulses  are  sometimes  sampled  incompletely,  an  alternative  to  matched 
hltering  is  desired.  The  threshold  crossing  technique  described  earlier  and  illustrated 
in  Fig.  5.11  is  an  acceptable  alternative  to  the  matched  hlter  for  this  LADAR  data 
set.  It  is  simple  and  performs  well  at  range  estimation. 

The  hrst  range  estimation  result  is  shown  in  Fig.  5.12a.  This  data  was  collected 
when  the  panel  was  tilted  32  degrees  from  the  LADAR  LOS.  Examples  of  the  raw  and 
equalized  data  used  in  these  range  calculations  are  shown  in  Figs.  5.8b  and  5.10b, 
respectively. 

The  raw  range  data  in  Fig.  5.12a  does  not  match  the  true  tilt  of  the  target  panel. 
The  155  cm  high  panel,  when  tilted  32  degrees,  has  a  depth  of  82  cm.  However,  the 
gain  variation  makes  the  target  depth,  as  measured  before  equalization  (shown  in 
blue),  appear  to  be  about  1.2  m.  The  range  estimates  made  with  the  equalized  data 
(shown  in  red  in  Fig.  5.12a)  do  not  match  the  target  panel’s  tilt  exactly.  However, 
they  are  much  more  precise  than  the  raw  range  measurements.  The  root  mean  squared 
error  (RMSE)  for  the  gain  equalized  range  estimates  was  8  cm.  The  RMSE  for  the 
raw  data  was  13  cm. 

In  Fig.  5.12b,  the  target  panel  was  tilted  22  degrees  from  LADAR  LOS.  In  this 
case,  gain  variation  in  the  raw  data  still  causes  a  dramatic  overestimate  of  the  panel’s 
tilt.  The  RMSE  is  15  cm.  The  equalized  data  matched  the  slope  of  the  panel  almost 
exactly.  The  gain  equalized  range  measurements  have  an  RMSE  of  only  4  cm.  This 
error  is  mostly  due  to  noise  in  the  laser  signals  themselves  rather  than  gain  variation 
distortions. 

In  Fig.  5.12c,  the  target  panel  was  oriented  perpendicular  to  the  LADAR’s 
line-of-sight.  The  range  error  for  the  raw  data  was  10  cm.  However,  when  the  gain 
equalized  data  was  processed,  the  RMSE  dropped  to  4  cm.  Because  the  target  is 
at  basically  the  same  range  in  all  target  pixels,  the  laser  pulse  distortion  due  to 
gain  variation  should  be  roughly  equal  on  all  signals.  However,  there  is  still  more 
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(c)  (d) 

Figure  5.12:  (a)  Target  panel  tilt  is  32  degrees,  (b)  Target  panel  tilt  is  22  degrees, 

(c)  Target  panel  tilt  is  0  degrees,  (d)  Target  panel  tilt  is  -20  degrees. 
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Table  5.1:  Range  Error  in  APD  Ar¬ 
ray  LADAR  Data. 


Panel  Tilt 

Raw  Data  RMSE 

Equalized  Data  RMSE 

-20  deg. 

9  cm 

6  cm 

0  deg. 

10  cm 

4  cm 

9  deg. 

14  cm 

3  cm 

15  deg. 

15  cm 

2  cm 

22  deg. 

15  cm 

4  cm 

27  deg. 

14  cm 

6  cm 

32  deg. 

13  cm 

8  cm 

fluctuation  in  the  range  estimates  made  with  the  raw  data  than  with  the  equalized 
data. 

The  hnal  example  of  range  estimation  is  illustrated  in  Fig.  5.12d.  In  this  case 
the  target  panel  is  -20  degrees  from  LOS.  Raw  and  gain  equalized  data  for  this  test 
are  illustrated  in  Figs.  5.8a  and  5.10a.  The  sign  change  of  the  panel  tilt  causes  the 
range  errors  to  occur  opposite  of  how  they  appear  in  Figs.  5.12a-c  where  the  tilt 
angle  was  nonnegative.  However,  the  gain  variation  technique  still  helped.  For  range 
estimates  made  with  raw  data,  the  RMSE  is  9  cm.  Using  gain  equalization,  the  range 
error  drops  to  6  cm. 

The  range  RMSEs  for  the  entire  data  collection  are  listed  in  Table  5.1.  For 
target  panel  tilt  angles  between  0  and  27  degrees,  the  raw  data  RMSE  is  at  least 
double  the  RMSE  for  range  measurements  made  with  the  equalized  data.  For  the 
large  magnitude  tilt  angles  of  -20  and  32  degrees,  the  raw  data  range  RMSE  is  higher 
than  for  the  equalized  data,  but  only  by  about  50%.  Thus,  there  appears  to  be  a  limit 
on  the  amount  target  tilt  that  the  gain  equalization  algorithm  can  tolerate  before  its 
performance  degrades  signihcantly. 
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5-4  Chapter  Summary 

Avalanche  photodiodes  are  used  in  LADAR  systems  to  detect  laser  light  reflected 
from  targets.  Arrays  of  APDs  can  be  very  dense  and  compact,  but  can  be  challenging 
to  work  with  because  of  variations  in  the  relationship  between  their  input  and  output 
signals.  In  Sec.  5.1,  LADAR  images  collected  by  a  system  using  an  APD  array  as  a 
photodetector  were  shown.  The  effect  of  gain  variation  on  the  recorded  signals  was 
illustrated  in  that  section. 

In  Secs.  5.2  and  5.3,  methods  of  gain  equalization  were  proposed  and  tested.  In 
Sec.  5.2  a  multidimensional  Gauss-Newton  method  was  used  successfully  on  simulated 
LADAR  data.  However,  this  technique  did  not  work  well  on  real  LADAR  data.  In 
Sec.  5.3,  an  equalization  method  that  used  target-free  pixels  was  successfully  demon¬ 
strated.  Through  the  use  of  a  target  panel  in  LADAR  images,  it  was  shown  that 
gain  variation  could  introduce  significant  bias  into  range  measurements.  The  gain 
equalization  method  tested  in  Sec.  5.3  decreased,  but  did  not  eliminate,  the  range 
bias  in  the  three-dimensional  LADAR  images. 
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VI.  Conclusion 


All  of  the  research  in  this  dissertation  applies  to  the  problem  of  making  measure¬ 
ments  of  range  to  a  target  using  a  LADAR  system.  However,  range  measurements 
are  affected  by  many  different  factors.  Noise  models,  signal  shapes,  and  detector  elec¬ 
tronics  all  influence  the  precision  of  LADAR  range  estimates.  Thus,  a  wide  variety  of 
subjects  were  discussed  in  the  dissertation. 

6.1  Poisson-Distributed  LADAR  Signal  Models 

In  Chapter  III,  three  different  LADAR  signal  models  were  examined.  These 
models  are  all  variations  of  work  previously  published  in  [4].  That  paper  assumes  that 
laser  pulses  received  by  LADAR  systems  can  be  modeled  as  a  truncated  parabola  plus 
a  background  noise  level.  It  is  also  assumed  that  the  LADAR  system  is  shot-noise 
limited,  which  means  that  the  signal  and  noise  follow  the  Poisson  distribution. 

In  Sec.  3.1,  the  laser  pulse  model  from  [4]  was  revisited.  The  Cramer- Rao  lower 
bounds  for  range  and  other  signal  parameters  derived  in  that  paper  were  expressed 
in  terms  of  summations.  In  this  dissertation,  analytic  approximations  to  those  sum¬ 
mations  were  derived. 

LADAR  systems  are  sometimes  used  to  image  targets  beneath  tree  canopies.  In 
Sec.  3.2,  the  signal  model  is  generalized  to  include  two  signal  pulses.  The  first  signal 
represents  the  tree  canopy  and  the  second  is  the  target  on  the  ground.  CRLBs  were 
derived  for  this  model.  It  was  shown  that  the  presence  of  obscuration  did  not  affect 
the  range  CRLB  for  the  obscured  target. 

Some  targets  observed  by  LADAR  systems  are  extended.  That  is,  their  depth 
is  significant  compared  to  the  length  of  the  transmitted  laser  pulse.  In  this  case,  the 
reflection  process  varies  the  length  of  the  received  pulse.  It  had  been  assumed  that  the 
pulse  width  in  the  received  signal  was  known.  But  in  Sec.  3.3,  the  width  was  treated 
as  an  unknown  for  the  purpose  of  CRLB  derivations  in  order  to  model  variations  in 
target  depth.  It  was  shown  that  the  range  CRLB  is  the  same  regardless  of  whether 
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the  pulse  width  was  known  or  unknown.  However,  the  CRLB  is  still  dependent  on 
the  pulse  width  in  this  case. 

It  was  shown  in  Chapter  III  that  the  Poisson  noise  model  used  to  simulated 
shot-noise  limited  signals  made  calculation  of  maximum  likelihood  estimates  difficult. 
However,  in  most  cases  it  was  straightforward  to  create  a  relatively  simple,  though 
suboptimal  estimator  which  nearly  achieved  the  CRLB.  This  fact  was  verihed  through 
the  use  of  Monte  Carlo  simulations. 

6.2  Range  Precision  of  LADAR  Systems 

Chapter  IV  of  this  dissertation  focused  on  the  topic  of  the  precision  of  LADAR 
range  measurements.  Some  of  the  work  in  this  chapter  used  the  Poisson-distributed 
signal  models  developed  in  Chapter  HI.  Other  analysis  relied  on  signal  models  that 
assumed  additive  white  Gaussian,  rather  than  Poisson  noise.  The  use  of  laser  pulse 
shapes  other  than  the  truncated  parabola  was  also  discussed. 

In  Sec.  4.1,  the  range  CRLB  derived  in  Chapter  HI  was  discussed.  The  range 
CRLB  for  a  truncated  parabolic  pulse  mixed  with  Gaussian  rather  than  Poisson  noise 
was  derived.  It  was  shown  that  the  Gaussian-noise  CRLB  is  always  less  than  the 
Poisson-noise  CRLB.  Thus,  range  measurements  are  more  precise  in  the  presence  of 
Gaussian  noise  than  when  the  signal  is  mixed  with  shot  (Poisson)  noise. 

The  topic  of  pulse  averaging  was  also  discussed  in  Sec.  4.1.  It  was  also  shown 
that  for  a  system  that  can  divide  a  hxed  amount  of  laser  energy  into  multiple  pulses, 
that  that  range  precision  of  a  shot-noise  limited  system  is  maximized  when  all  energy 
is  transmitted  in  a  single  pulse.  That  is,  averaging  decreases  range  precision  unless 
the  total  amount  of  pulse  energy  is  increased. 

In  Sec.  4.2,  it  was  noted  that  some  LADAR  systems  operate  at  ranges  of  several 
kilometers.  Airborne  systems  in  particular  tend  to  operate  at  long  ranges.  Airborne 
LADAR  system  often  observe  a  scene  from  an  oblique  angle.  Unless  a  target  is 
oriented  normal  to  the  LADAR  line-of-sight,  the  laser  pulse  reflection  process  will 
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distort  the  shape  of  the  received  laser  signal.  It  was  shown  that  if  a  nanosecond- class 
dnration  beam  that  is  a  few  inches  in  diameter  at  the  transmitter  propagates  more 
than  a  few  kilometers,  the  effect  of  target  orientation  on  range  precision  can  not  be 
neglected. 

In  Sec.  4.2,  the  shape  of  the  received  laser  signal  was  derived  as  a  function  of 
target  snrface  orientation  and  other  system  parameters.  The  CRLB  was  calculated 
as  a  function  of  target  tilt.  The  effect  of  mitigation  strategies  such  as  increasing 
the  transmitted  beam  size  and  decreasing  the  pulse  duration  was  studied  and  the 
limitations  of  those  method  were  quantified.  This  analysis  assumed  a  laser  pulse  with 
a  Gaussian  shape  and  a  received  signal  mixed  with  AWGN. 

6.3  Gain  Variation  in  APD  Arrays 

The  topic  of  gain  variation  was  studied  in  Ghapter  VI.  Gompact  avalanche 
photodiode  arrays  are  being  used  to  create  three-dimensional  LADAR  images  with 
large  numbers  of  pixels.  Ideally,  the  amount  of  amplification,  or  gain,  in  an  APD 
would  be  constant.  However,  it  has  been  shown  in  this  dissertation  and  elsewhere 
than  gain  in  APD  arrays  varies  over  time.  Gain  variation  distorts  the  recorded  laser 
signals  and  introduces  bias  into  measurements  of  target  range. 

In  the  course  of  this  dissertation  research,  a  LADAR  data  collection  experiment 
was  done  to  examine  the  effect  of  gain  variation  (see  Sec.  5.1).  Three-dimensional 
images  were  collected  using  an  Air  Force  Research  Laboratory  LADAR  system.  The 
images  were  of  a  large  target  panel.  The  panel  was  tilted  at  various  orientations  with 
respect  to  the  LADAR  line-of-sight.  Examination  of  the  effect  of  target  tilt  clearly 
illustrated  the  effect  of  gain  variation  on  range-resolved  LADAR  data.  It  was  shown 
that  gain  variation  introduced  bias  into  range  measurements. 

Gain  variation  degraded  the  ability  of  LADAR  systems  to  make  range  measure¬ 
ments.  Therefore,  a  method  of  gain  equalization  that  corrects  for  variations  in  the 
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array  of  APDs  was  desired.  In  this  dissertation,  two  methods  of  gain  equalization 
were  proposed  and  tested. 

The  hrst  method  of  gain  equalization  was  developed  in  Sec.  5.2.  This  method 
was  based  on  the  Gauss-Newton  method  of  parameter  estimation.  A  model  was 
developed  for  range-resolved  LADAR  signals  observed  by  an  array  of  avalanche  pho¬ 
todiodes  in  the  presence  of  gain  variation.  This  model  assumed  Gaussian-shaped  laser 
pulses  and  white  Gaussian  noise.  The  method  was  tested  successfully  on  simulated 
data.  However,  the  Gauss-Newton  method  did  not  work  well  with  real  LADAR  data 
collected  from  the  target  panel.  This  relatively  sophisticated  approach  relied  on  too 
many  assumptions,  such  as  exact  knowledge  of  the  shape  of  the  signals  embedded  in 
the  noisy  data. 

The  second  method  of  gain  equalization  was  derived  and  tested  in  Sec.  5.3. 
This  method  used  averaged  signals  from  target-free  pixels  to  measure  the  change  in 
gain  over  time  in  the  APD  array.  This  background-averaging  method  was  tested  using 
data  collected  from  the  target  panel  tilt  tests.  It  was  shown  that  the  gain  equalization 
method  based  on  background  data  reduced,  though  did  not  eliminate,  the  bias  in  the 
range  measurements  to  the  target  panel.  Though  this  method  did  not  eliminate  the 
problems  introduced  by  gain  variation  in  APD  arrays,  it  had  the  advantage  of  being 
quite  simple.  This  method  could  potentially  be  implemented  in  near  real-time  in  a 
LADAR  receiver’s  ROIG. 

6.4  Summary 

Measurement  of  range  to  a  target  is  one  of  the  most  fundamental  problems 
in  the  held  of  remote  sensing.  Short  duration  laser  pulses  make  LADAR  systems 
ideal  instruments  for  precise  range  measurement.  This  dissertation  discussed  limits 
on  range  precision  of  LADAR  measurements  due  to  stochastic  noise  sources,  the  laser 
pulse  rehection  process,  and  characteristics  of  current  photodetectors.  The  work  can 
be  used  to  model  the  performance  of  current  LADARs  and  could  help  design  future 
LADAR  systems. 
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Appendix  A.  Photoelectron  Statistics 

Much  of  the  analysis  in  this  dissertation  relies  on  the  Poisson  approximation.  That  is, 
it  is  frequently  assumed  that  probability  distribution  describing  number  of  photoelec¬ 
trons  observed  in  an  imaging  LADAR  has  a  Poisson  distribution.  In  this  appendix, 
the  origin  of  and  justification  for  this  assumption  is  discussed.  Some  of  the  derivations 
used  in  the  Poisson  approximation  are  reviewed. 


A.l  Integration  of  Intensity 

One  of  the  early  references  on  photon  statistics  is  ^Some  Effects  of  Target- 
Induced  Scintillation  on  Optical  Radar  Performance”  by  Goodman  [11].  That  paper 
notes  that  the  total  energy  W  incident  on  an  optical  system’s  receiving  aperture  is 


W  = 


Ar(x,  y)w(x,  y)dxdy 


—  OC  —OO 


(A.l) 


where  x  and  y  are  the  aperture  coordinates,  w{x,y)  is  the  energy  density  in  the 
aperture  plane,  and  Affx^y)  is  the  binary  receiver  aperture  function.  For  a  circular 
aperture  with  diameter  Da,  this  function  is 


Ar{x,y) 


1  when  <  Da/2 

0  when  ySphfffff  >  Da/2 


(A.2) 


Using  results  from  earlier  works  [33,  34]  Goodman  notes  in  Eq.  15  of  [11]  that  the 
total  energy  is  approximately  a  gamma  random  variable.  Therefore,  the  PDF  of  W 

is 


p{W) 


1 

r(A^) 


aMwM-l^-aW 


when 


lU  >  0 


0  when  lU  <  0 


(A.3) 


where  a  and  A4  are  related  so  that  the  mean  energy,  which  is  denoted  W,  is 


W  = 


M 

a 


(A.4) 
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The  variance  is  the  energy  is 


Goodman  shows  that 


r  r  r)  • 

(A.5) 

M 

(A.6) 

a  =  ^ 

W 

and  that  M  is  the  diversity  or  the  nnmber  of  spatial  correlation  cells  snbtended  by 
the  apertnre. 


The  second  moment  of  the  energy  can  be  evalnated  by  integrating  the  antocor- 
relation  of  the  energy  density  fnnction 


1T2  = 


Ar(xi,yi)  Ar  (X2 , 2/2)  (^1 , 2/i ;  Xa ,  2/2 )  dxidpi  dx2dy2  •  (A.  7) 


—  OO  — OO  — OO  — OO 


This  integral  reqnires  the  energy  density  antocorrelation  //i;  X2, 2/2),  which  is 

(Eq.  12  of  [11]) 

RU^i,yi;x2,y2)  =  E[w{xi,yi)w{x2,y2)]  ■  (A.8) 

The  diversity  parameter  is  (Eq.  18  of  [11]) 

M  =  ^  ^ -  (AJ) 

f  f  f  f  Ar(xi,yi)Ar(x2,y2)j7(^i,yi',^2,y2)j^dxidyidx2dy2 

—  OO  — OO  — OO  — OO 

where  7(xi,  2/1;  2:2, 2/2)  is  the  complex  degree  of  coherence  (see  Eq.  5.2-11  of  [12]). 


A. 2  Photon  Arrival  Statistics 

It  is  assnmed  that  arrival  of  photons  is  a  Poisson  process  [25].  Given  incident 
energy  W,  the  distribntion  of  the  nnmber  of  photons  arriving  in  some  time  interval 
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is  a  Poisson  random  variable. 


p{k\W)  = 


where 


(A.IO) 


«:  =  0,1,2,.... 


where  k  is  the  nnmber  of  photons  arriving  and  N  is  the  average  nnmber  of  photons 
observed.  The  mean  and  the  energy  W  are  related  by 


(A.ll) 


In  the  above  eqnation,  h  is  Planck’s  constant  and  v  is  the  optical  freqnency. 

The  distribntion  of  photons  arriving  conditioned  on  the  energy  is  given  by  Eq. 
A.IO.  However,  the  energy  W  is  itself  an  random  variable.  The  joint  distribntion  of 
k  and  W  is  jnst  the  prodnct  of  the  energy  density  from  Eq.  A. 3  and  the  conditional 
density  from  Eq.  A.IO.  The  joint  density  is 


p{k,W)  =  p{k\W)P{W) 


(A.12) 


The  marginal  density  of  k  is  obtained  by  integrating  over  joint  density  p{k,  W). 


OO 


(A.13) 


0 


It  was  shown  in  Eqs.  28-30  of  [11]  that  this  is 


where 


(A.14) 


A;  =  0, 1,  2, ... . 
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Negative  Binomial  and  Poisson  PDFs  Negative  Binomial  and  Poisson  PDFs 


k  k 


(a)  (b) 


Figure  A.l:  (a)  Negative  binomial  densities  with  M.  =  1  and  Poisson  density, 
(b)  Negative  binomial  densities  with  M.  =  2  and  Poisson  density. 


This  is  a  negative  binomial  density  with  mean  N  and  diversity  parameter  M..  The 
diversity  is  given  in  Eq.  A. 9. 


A. 3  Poisson  Approximation  to  Negative  Binomial  Distribution 

The  mean  number  of  photons  per  speckle  correlation  cell  is  called  the  count 
degeneracy  parameter  [13].  It  is  denoted  6c  ■ 


6c  = 


N 

M' 


(A.15) 


It  was  noted  in  [11]  that  when  the  count  degeneracy  parameter  becomes  small 


6c  «  1 


(A.16) 


that  the  negative  binomial  density  converges  to  the  Poisson  distribution  with  mean 
N.  This  characteristic  is  illustrated  in  the  Fig.  A.l.  The  plots  in  Fig.  A.l  show 
that  as  the  mean  number  of  photons  per  speckle  cell  approaches  zero,  the  negative 
binomial  density  approaches  the  Poisson  density. 


160 


0.1 

Q. 

0.05 

(a)  (b) 

Figure  A. 2:  (a)  Negative  binomial  densities  with  5c  =  0.5  and  Poisson  densities, 

(b)  Negative  binomial  densities  with  5c  =  0.75  and  Poisson  densities. 

The  negative  binomial  density  also  converges  to  the  Poisson  density  as  the 
diversity  goes  to  inhnity  [5].  In  this  case  the  convergence  is  very  slow  and  is  difficult 
to  illustrate.  In  Fig.  A. 2,  the  ratio  of  the  mean  number  of  photons  to  the  diversity 
5c  =  N / M.  is  not  near  zero.  Therefore,  the  condition  required  by  Eq.  A.  16  is  not 
satished.  However,  as  M.  gets  large,  the  negative  binomial  density  is  slowly  converging 
to  the  Poisson  density. 

A. 4-  Photon  Arrivals  and  Quantum  Efficiency 

LADAR  systems  rely  on  photodetectors  to  convert  photons  to  a  measurable 
signal.  The  process  of  converting  photons  to  electrons  (typically  called  photoelectrons) 
is  imperfect.  That  is,  some  photons  incident  on  the  detector  will  fail  to  convert  to 
signal.  The  fraction  of  photons  converted  to  photoelectrons  is  called  the  quantum 
efficiency  and  is  denoted  y. 

The  statistics  describing  photon  arrival  densities  were  derived  using  mathemat¬ 
ics  that  model  optical  propagation.  Those  derivations  tell  how  many  photons  are 
incident  on  a  detector  in  some  time  interval.  However,  only  some  of  the  incident 
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photons  are  successfully  converted  to  photoelectrons  by  that  detector.  Undetected 
photons  are  not  of  interest.  It  is  the  photoelectron  density  that  is  desired. 

In  ‘^Effect  of  random  deletion  and  additive  noise  on  bunched  and  anitbunched 
photon- counting  statistics” ,  Teich  and  Saleh  [50]  considered  several  different  photon 
arrival  densities.  They  assumed  that  the  probability  of  each  photon  being  converted 
to  a  photoelectron  was  hxed  and  statistically  independent  of  the  conversion  process  for 
all  other  arriving  photons.  They  proved  that  if  photon  arrivals  are  Poisson  distributed 
with  mean  N,  that  the  photoelectron  density  is  pN.  For  a  negative  binomial  density 
with  mean  N  and  diversity  Ai,  the  photoelectron  density  is  also  negative  binomial 
with  mean  rjN.  They  also  showed  that  the  diversity  of  the  photoelectron  density  is 
the  same  as  the  diversity  of  the  photon  density. 
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Appendix  B.  Parabolic  and  Gaussian  Pulse  Properties 


This  appendix  contains  a  discnssion  of  the  dnration  and  bandwidth  of  the  temporal 


LADAR  signal  models  nsed  in  this  dissertation.  The  parabolic  pulse  shapes  that  are 
used  in  Chapters  II  and  III  are  discussed.  The  Gaussian  pulse  models  from  Chapters 
IV  and  V  are  also  reviewed.  The  full-width  half-maxima  of  these  pulses  are  calculated. 
The  energy  spectral  densities  are  derived  using  the  Fourier  transforms  of  the  pulses. 


The  durations  and  bandwidths  of  parabolic  and  Gaussian  pulses  are  compared. 


The  signals  described  in  this  appendix  are  the  photodetector  output  signals. 
The  energy  or  power  in  these  signals  is  generally  not  equal  to  the  optical  energy  or 
power  in  the  received  laser  pulse.  The  photocurrent  signals  are  measured  in  units  of 
normalized  current  or  normalized  voltage.  These  normalized  units  are  proportional 
to  the  received  optical  power  observed  by  the  LADAR  photodetector.  These  relation¬ 
ships  are  explained  in  more  detail  in  Sec.  B.l. 

B.l  Parabolic  Pulse  Properties 

The  parabolic  pulse  equation  centered  at  time  zero  is 


(B.l) 


The  parabolic  pulse’s  FWHM  is  given  by  Eq.  4.45.  The  Fourier  transform  of  the 
parabolic  pulse  is 


Pw 


(B.2) 


G 


Ap^G 


sin  {2np^f)  -  {27ip^f)  cos  Gnp^f) 
{27rp^ff 
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Parabolic  Pulse 


(a) 


(b) 


Figure  B.l;  (a)  Parabolic  pulse  with  width  =  10  ns  and  peak  normalized  voltage 
G  =  2  W^/2 

(b)  Pulse  Fourier  transforms  for  widths  =  5,  10,  and  20  ns. 


A  10  ns  parabolic  pulse  is  illustrated  in  Fig.  B.la.  Fourier  transforms  for  5,  10, 
and  20  ns  pulses  are  shown  in  Fig.  B.lb.  The  pulses’s  units  are  normalized  voltage, 
which  is  the  voltage  divided  by  the  square  root  of  the  resistance  (e.g.  Volts  per  root 
Ohm).  The  units  of  normalized  voltage  are  equal  to  the  square  root  of  the  power  (e.g. 
Watts). 


VPower  = 


The  reason  for  using  normalized  voltage  in  place  of  regular  voltage  is  that  normalized 
voltage  produces  energy  and  power  spectra  with  the  correct  units  of  energy  times 
time  and  power  times  time,  respectively. 

Alternately,  the  photocurrent  could  be  analyzed  in  terms  of  normalized  current. 
Normalized  current  is  current  times  the  square  root  of  resistance  (e.g.  Amps  times 


^Voltage  X  Current 


Voltage 

Resistance 


Voltage  X 

Voltage 

■y/Resistance 


(B,3) 
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root  Ohms). 


VPower  =  ■\/Voltage  x  Current 


=  VCurrent  x  Current  x  Resistance 


=  Current  x  \/Resistance. 


(B.4) 


As  with  normalized  voltage,  normalized  current  will  produce  energy  and  power  spectra 
with  the  correct  units  of  energy  times  time  and  power  times  time. 

A  signal’s  energy  spectral  density  (ESD)  is  the  magnitude  squared  of  its  Fourier 
transform  [32].  The  ESD  of  the  parabolic  pulse  from  Eq.  B.l,  which  is  denoted  Gp{f), 

is 


-2  'sin  (27rp^/)  -  (27rj9^/)  cos  (27rp^/)'  ^ 
{^Trpwff 


Gpif)  =  \Ppif)\^  =  IGplG' 


(B.5) 


The  ESDs  for  parabolic  pulses  with  various  widths  are  illustrated  in  Fig.  B.2a. 

In  most  LADAR  systems,  the  received  signal  is  digitally  sampled.  In  a  sampled 
signal,  pulse  energy  inside  the  Nyquist  frequency,  which  is  half  the  sampling  frequency 
ifn  =  fs/‘P)i  is  sampled  properly.  However,  energy  outside  than  the  Nyquist  frequency 
is  aliased  [30].  Aliasing  will  cause  interference  and  will  degrade  the  performance  of 
parameter  estimation  techniques.  The  fraction  of  pulse  energy  that  is  aliased  can  be 
calculated  by  integrating  the  ESD  over  frequencies  with  magnitudes  greater  than  /„, 
then  dividing  by  the  total  energy.  Since  the  ESD  has  even  symmetry,  the  fraction  of 
energy  outside  the  Nyquist  frequency  is 


J  oo  f„ 

5  G,(f)df  f  a,{f)df  -  f  G/f)df 


fn  0  0 


(B.6) 


oo  oo 


0 


0  0 
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where  the  energy  in  the  parabolic  pulse  Ep  is 


/  G,(f)df=  /  \P,(f)fdf. 


(B.7) 


An  efficient  way  to  calculate  the  pulse  energy  is  through  the  use  of  Parseval’s  theorem, 
which  states  that 


(B.8) 


Combining  Eqs.  B.l,  B.7,  and  B.8  yields 


Pw 


Ep  = 


Gll-^ 

pI 


16 

dt  =  —G'^p^ 
15 


(B.9) 


-Pw 


Using  this  result  in  Eq.  B.6  shows  that  the  fraction  of  energy  that  is  greater  than  the 
Nyquist  frequency  can  be  written 


fn 


Apifn)  =  1  -  30p^ 


sin  (27rp^/)  -  (27rp^/)  cos  {27ip^f) 
{‘27rpwf)^ 


df.  (B.IO) 


The  fraction  of  a  parabolic  pulse’s  energy  outside  the  Nyquist  frequency  is  il¬ 
lustrated  in  Fig.  B.2b.  For  the  20  ns  pulse,  about  3  x  10“^%  of  the  pulse  energy  is 
aliased  when  the  sampling  rate  is  one  Gigahertz  (/„  =  500  MHz).  However,  for  the  5 
ns  pulse,  about  2  x  10“^%  of  the  energy  is  aliased. 


B.2  Gaussian  Pulse  Properties 

Let  the  Gaussian  pulse  equation  is 

Pg{t)  =  Gexp(^-^^  (B.ll) 
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ESD  of  Parabolic  Pulse 


Frequency  (GHz) 


(a) 


(b) 


Figure  B.2;  (a)  Parabolic  pulse  energy  spectrum  for  pulse  widths  =  5,  10,  and 

20  ns. 

(b)  Fraction  of  pulse  energy  outside  Nyquist  frequency. 


where  G  is  the  gain  and  w  is  the  variable  describing  the  width  of  the  pulse.  At 
distance  w  from  the  pulse  center,  the  Gaussian  pulses’s  amplitude  is  times  its 
peak  height.  The  full-width  half- maximum  of  this  pulse  is  given  by  Eq.  4.36. 

The  Fourier  transform  of  the  parabolic  pulse  is 

OO 

—  OO 

The  Gaussian  pulse  and  its  Fourier  transform  are  illustrated  in  Fig.  B.3.  The  ESD 
of  the  Gaussian  pulse,  which  is  denoted  Gg{f),  is 


Ggif)  =  \PM)f  =  (B.13) 


The  energy  in  the  Gaussian  pulse,  denoted  Eg,  is 


j  Gg{f)dt  =  ^wG^. 


—  OO 


(B.14) 
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The  fraction  of  energy  that  is  outside  the  Nyquist  frequency,  denoted  Ag{fn),  is 


Mfn) 


f  G,{f)df 


fn 


I  G,i  f  )df 

0 

fn 

1  -  /  G,{f)df 


1  - 


^-(27Twfy 


df 


2-KWfn 

1  -  ^ 


1  -  erf(2irw/„) 


(B.16) 


where  the  change  of  variable  u  =  2'Kwf  was  used  and 


X 


(B.16) 


is  the  error  function  (see  Eq.  7.1.1  of  [1]).  The  ESD  of  the  Gaussian  pulse  and  energy 
fraction  are  illustrated  in  Fig.  B.4.  The  ESD  of  the  Gaussian  pulse  decays  much 
faster  than  the  ESD  for  the  parabolic  pulse  (shown  in  Fig.  B.2a).  Similarly,  the 
function  describing  the  fraction  of  Gaussian  pulse  energy  Ag{fn)  decays  very  quickly 
compared  to  the  same  function  for  the  Gaussian  pulse. 


B.3  Parabolic  and  Gaussian  Pulse  Comparison 

The  parabolic  and  Gaussian  pulses  are  compared  in  Fig.  B.5.  Both  pulses  have 
a  10  ns  FWHM  and  1/iJ  of  energy.  The  gain  for  the  pulses  and  the  FWHM  are 
almost  identical.  However,  the  parabolic  pulse  is  narrower  (see  Fig.  B.5a)  and  its 
energy  fraction  that  falls  outside  a  Nyquist  frequency  is  much  higher  (see  Fig.  B.5b). 
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Figure  B.3:  (a)  Gaussian  pulse  with  width  w  =  5  ns  and  peak  normalized  voltage 

G  =  2  W^/2 

(b)  Pulse  Fourier  transforms  for  widths  w  =  2,  5,  and  10  ns. 


ESD  of  Gaussian  Pulse 


(a) 


Figure  B.4:  (a)  Gaussian  pulse  energy  spectrum  for  pulse  widths  w  =  2,  5,  and  10 

ns. 

(b)  Fraction  of  pulse  energy  outside  Nyquist  frequency. 
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Figure  B.5:  (a)  Parabolic  and  Gaussian  pulses  with  10  ns  FWHM  and  1/iJ  energy, 

(b)  Fraction  of  pulse  energy  outside  Nyquist  frequency. 


Power  Full- Width  Half-Maximum 


The  FWHM  can  be  measured  from  the  pulse’s  power  instead  of  the  pulse’s 
normalized  voltage.  The  power  is  the  normalized  voltage  squared.  The  power  FWHM 
for  a  parabolic  pulse  is 


6p 


(B.17) 


The  power  FWHM  for  a  Gaussian  pulse  is 


(5j  =  2w  [log(2)|'''^ . 


(B.18) 


The  power  FWHMs  for  parabolic  and  Gaussian  pulses  are  illustrated  in  Fig.  B.6. 


B.5  Time- Bandwidth  Products 

The  time-bandwidth  product  (TBP)  of  a  waveform  is  used  as  a  measure  of  its 
quality.  It  is  a  measure  of  the  compactness  of  a  waveform  in  both  time  and  frequency. 
The  TBP  is  computed  in  the  same  way  as  beam  quality  was  calculated  in  Ghapter  H. 
It  is  the  product  of  the  mean  square  duration  at  and  the  mean  square  bandwidth  cx/ 
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Parabolic  Pulses 


Gaussian  Pulses 


Figure  B.6:  (a)  Power  FWHM  for  parabolic  pulses, 

(b)  Power  FWHM  for  Gaussian  pulses. 


of  a  waveform  [52]. 

TBP  =  ataj.  (B.19) 

Given  arbitrary  waveform  p{t),  the  quantities  needed  to  calculate  the  time-bandwidth 
product  are 

OO 

/  t^p^{t)dt 

2  -OO 

CXD 

/  p^{t)dt 

—  OO 

and 

OO 

I 

9  — OO 

V  =  - 

I  IP(f)l  df 

—  OO 

where  P(f)  is  the  Fourier  transform  of  the  waveform. 

The  waveform  uncertainty  relation  is  a  lower  bonnd  on  the  time-bandwidth 
prodnct  of  any  waveform  (see  p.  600  of  [52]).  The  waveform  nncertainty  relation  is 

atCTf  >  (B.22) 


(B.20) 


(B.21) 
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For  the  Gaussian  pulse  from  Eq.  B.ll,  the  mean  square  duration  squared  is 


w 


a 


t—gaussian 


and  the  mean  square  bandwidth  squared  is 


a 


f-gaussian  g-jj-2yj2 ' 


(B.23) 


(B.24) 


Thus,  the  time-bandwidth  product  for  a  Gaussian  pulse  is 


TBP 


gaussian 


^  0.0796. 

4:71 


(B.25) 


Gomparison  of  this  result  with  the  waveform  uncertainty  relation  (Eq.  B.22)  shows 
that  the  Gaussian  waveform  achieves  the  bound.  It  is  the  only  waveform  that  achieves 
this  bound. 

Given  the  parabolic  pulse  from  Eq.  B.l,  the  mean  square  duration  squared  is 


a 


2 

t— parabolic 


7 


and  the  mean  square  bandwidth  squared  is 


a 


2 

f— parabolic 


5 


(B.26) 


(B.27) 


Thus,  the  time-bandwidth  product  for  a  parabolic  pulse  is 


TBP 


parabolic 


1 


14  J 


0.0951. 


(B.28) 


Gomparison  of  the  time-bandwidth  products  of  the  Gaussian  pulse  (Eq.  B.25) 
and  the  parabolic  pulse  (Eq.  B.28)  shows  that,  as  expected,  the  parabolic  pulse’s 
time-bandwidth  product  is  greater.  However,  the  difference  is  not  that  great.  The 
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ratio  is 


1/2 


TBP 


parabolic 


gaussian 


^  1.195. 


(B.29) 


Thus,  the  parabolic  pulse’s  TBP  is  less  than  20%  higher  than  for  the  Gaussian  pulse. 
Therefore,  the  truncated  parabolic  pulse  is  a  relatively  compact  waveform. 


The  mean  squared  duration  of  the  Gaussian  pulse  in  terms  of  the  FWHM  is 
(see  Eqs.  4.36  and  B.23) 


(Jt-gaussian  =  — - ^—77^  ~  0.3007^.  (B.30) 

4[log(2)]'/2 

For  the  parabolic  pulse,  the  mean  squared  duration  is  (see  Eqs.  4.45  and  B.26) 

(^t-paraboiic  =  ~  0.267rp.  (B.31) 

The  duration  of  the  parabolic  pulse  is  shorter  than  the  duration  of  the  Gaussian  pulse 
by  over  10%.  The  parabolic  pulse  is  an  example  of  a  pulse  that  is  a  compact  in  the  time 
domain.  It  is  the  relatively  wide  bandwidth  of  the  parabolic  pulse,  which  is  illustrated 
in  Fig.  B.5,  that  makes  its  TBP  greater  than  the  Gaussian  pulse’s  time-bandwidth 
product. 
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Appendix  C.  Lower  Bound  on  Estimator  Varianee 

In  this  appendix,  it  is  demonstrated  that  the  reciprocal  of  any  diagonal  element  of  a 
Fisher  information  matrix  is  a  lower  bound  for  the  Cramer-Rao  lower  bound  associated 
with  that  parameter.  In  other  words,  the  FIM  diagonal  element’s  reciprocal  is  a  lower 
bound  on  the  variance  of  unbiased  parameter  estimates.  However,  it  is  generally  a 
looser  lower  bound  than  the  CRLB. 

The  fact  FIM  diagonal  element  reciprocals  are  lower  bounds  for  the  CRLB  has 
already  been  noted  by  Van  Trees  for  the  special  case  of  2x2  Fisher  information  matri¬ 
ces  [52].  However,  that  reference  does  not  discuss  cases  where  the  Fisher  information 
matrix  is  larger  than  2x2.  In  this  appendix  it  is  shown  that  the  property  is  true  for 
any  nonsingular  FIM,  regardless  of  its  size. 

C.l  Decomposition  of  the  Fisher  Information  Matrix 

Fisher  information  matrices  are  nonnegative  dehnite.  This  follows  from  the  fact 
that  they  can  be  dehned  as  covariance  matrices  [23].  It  is  possible  for  the  FIM  to  be 
singular.  In  that  case,  the  FIM  inverse  does  not  exist  and  there  are  no  unbiased  hnite 
variance  estimators  of  the  unknown  parameters  [49]. 

If  the  FIM  is  also  nonsingular,  then  it  is  positive  dehnite.  Any  positive  dehnite 
matrix  has  the  following  eigenvalue  decomposition  [36]. 

J  =  USU^  (C.l) 

where  U  is  a  unitary  matrix  and  S  is  diagonal.  The  notation  ^  denotes  Hermitian 

(conjugate)  transpose.  By  exploiting  the  fact  that  U  is  unitary  and  S  is  diagonal, 
Eq.  C.l  can  be  written 

M 

^  =  ^  (C.2) 

m=l 

where  M  is  number  of  rows  and  columns  in  the  FIM.  The  vector  is  the  m-th 
column  of  U. 

U  =  [ui  U2  ...  um]  •  (C.3) 
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The  diagonal  elements  of  S  are  denoted  A^-  Since  it  has  been  assnmed  that  the  FIM 
is  positive  definite,  all  of  its  eigenvalnes  are  positive  [18]. 

Am  >  0  V  m.  (C.4) 


The  inverse  of  the  Fisher  information  matrix  is 

M 

3  ^  — UmU^.  (C.5) 

m=l  ^ 

The  Cramer- Rao  lower  bonnd  for  estimates  of  the  n-th  parameter  is  the  reciprocal  of 
the  n-th  diagonal  element  of  the  inverse  of  the  FIM.  That  is, 


Var 


>  [J-i]  . 


Combining  Eqs.  C.5  and  C.6  prodnces 


Var 

On 

> 

M 


EJ- 

^UmU: 


jn=l 

M  . 


m=l 

M 


E  w  Pi„. 


m=l 


(C.6) 


(C.7) 


C.2  Relationship  Between  CRLB  and  FIM  Diagonal  Elements 

In  this  section,  the  relationship  between  the  Cramer-Rao  lower  bonnd  and  the 
diagonal  elements  of  the  Fisher  information  matrix  is  stndied.  Using  Eq.  C.2,  the 
n-th  diagonal  element  of  the  FIM  can  be  written 


M 


[Ji„„  =  E  [ui, 


m=l 


(C.8) 
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Let 


^n.rn. 


n.m  V  ^'TTi 


[U], 


and 


^n,m  ^ — 

V  '^m 


[u], 


Eqs.  C.8  and  C.9  combine  to  form 


M 


E 

m=l 


^*^]nn  ^  j  ^n,m' 


Eqs.  C.6,  C.7,  and  C.IO  combine  to  prodnce 


M 


[n‘]„„  =  5:<. 


m=l 


(C.9) 


(C.IO) 


(C.ll) 


(C.12) 


The  Canchy-Schwarz-Bnniakowsky  ineqnality  is  (p.  1049  of  [14]) 


(C.13) 


All  real  nnitary  matrices  are  orthogonal.  Becanse  U  is  orthogonal,  its  rows  and 
colnmns  are  each  orthonormal  bases  [23].  Using  this  fact,  along  with  Eqs.  C.9  and 
C.IO,  shows  that 


M  M 

'^n,mbn,m  =  |  [U]„ 

m=l  m=l 


1 


(C.14) 


Therefore, 


M 


> 


/  y  ^n,m  —  ^ 

m=l  fj2 

^n,r. 

m=l 


(C.15) 


When  Eq.  C.15  is  combined  with  Eqs.  C.ll  and  C.12,  the  following  ineqnality  is 
prodnced. 


(C.16) 
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Combining  this  result  with  the  CRLB  from  Eq.  C.7  produces  a  lower  bound  for  the 


variance  of  unbiased  estimates  of  parameter  On- 


Var 


(C.17) 


This  bound  is  generally  looser  than  the  CRLB.  However  it  is  of  value  when  there  is  a 
complicated  Fisher  information  matrix  that  is  cumbersome  to  invert.  In  that  case,  a 
simple  lower  bound  can  be  obtained  by  simply  taking  the  reciprocal  of  the  appropriate 
FIM  diagonal  element. 
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